What's new
==========

This page records user-visible changes to the **modern** RISE Toolbox --
new features, behaviour changes, and things you may need to know when
upgrading. Newest entries go on top.

.. note::

   This page covers the modern toolbox only. The legacy toolbox
   (``RISE_Tbx_Beta``) has its own release notes in ``rise-stable-docs``;
   the two are never mirrored into each other.

Unreleased
----------

- **Fixed: an observed driver of the transition probabilities in a
  model with occasionally-binding constraints.** Filtering such a model
  stopped with an index error: the constraint engine evaluated the
  transition matrix without the driver's value. It now receives it, so
  time-varying transition probabilities driven by data and constraints
  such as the zero lower bound can be filtered and estimated together.

- **Fixed: the first period of the constrained filter.** ``imm`` updated
  through the constrained step from inputs set at the end of each period
  for the next, and nothing set them for the first observed period: its
  likelihood ignored the constraints. They are now set from the initial
  distribution whenever it is the stationary one, so the first period is
  treated like any other. Models whose first period is far from their
  constraints are unchanged bit for bit; an exact diffuse start, or an
  initial state given as the first prediction, keep today's update.

- **Changed: the steady state of the lead of a growing variable moves with
  the trend.** In a model with a balanced growth path, ``stst(V{+1})``
  written in a dynamic equation was the value of the path at the expansion
  date, a constant: a hand-written anchor such as
  ``stst(C{+1})+mcd*(C{+1}-stst(C{+1}))`` then had no balanced growth path
  for ``mcd < 1`` and the solution drifted without shocks. It is now the
  moving reference of ``@cognitive_discount``: exogenous unit-root trends at
  their realized value, the other trends at their previous value grown at
  trend. For the log price level with ``p = p{-1} + pi``,
  ``stst(p{+1})+m*(p{+1}-stst(p{+1}))`` is ``p + pis + m*(pi{+1} - pis)``.
  Models without trends, stationary variables and ``stst()`` of current
  values or lags are unchanged. See :ref:`stst-of-growing-variables`.

- **Fixed: ``'bm'`` on the planner-discount knife-edge.** The solve failed
  (retcode 30) because the anchor condition has a continuum there; the
  solution is now the limit of those at nearby planner discounts. See
  :doc:`/ModelShapes/DSGE/Optimal policy`.

- **Fixed: the smoothed shocks of the first period.** They were exactly
  zero: the filter gave the shocks of the first period an initial covariance
  independent of the first prediction of the state, which they move. They
  are now the exact ones. On fs2000, RISE's smoothed shocks match Dynare's in
  every period; they differed by 6e-4 in the first. The likelihood and the
  filtered and smoothed variables do not change.

- **Changed: the known-regimes filter (``myKnownRegimesFilter``) takes
  the regime of each period, and refuses what it cannot filter.** In
  ``filter(m, kf_user_algo = {@myKnownRegimesFilter, regs})``, ``regs(t)``
  is the regime of period t: the regime of the solution that produces
  period t, whose standard deviations the shocks of period t have. The
  filter used to apply each regime one period late; a path shifted to
  compensate must now be passed as is. The constant of each period is
  that of its regime (models whose regimes have different steady states
  were filtered with the constants of the wrong regimes). The initial
  state and covariance are those of the regime of the first period
  instead of an ergodic average; ``kf_user_init`` overrides them.
  Constraints (unless ``simul_disable_constraints``), news shocks,
  ``solve_order > 1`` and deterministic exogenous variables raise errors
  with identifiers ``rise:knownRegimes:*``. With shock standard deviations
  written as switching parameters the filter is the heteroskedastic
  filter: see :ref:`known-regime-paths`.

- **Changed: the auxiliary equations of an ``@outside_optimization``
  block are exempt with it.** The equation RISE adds for a lead or lag
  beyond one inside a block (``LEAD_1_B{t} = B{t+1}`` for a ``B{t+2}``)
  was a constraint of every player, with a multiplier zero along the
  solution, and could raise the feedback warning. It now belongs to the
  block, in both forms: the exempt players do not include it in their
  system and take its variable as given (no derivative with respect to
  it). The solutions are unchanged; the multipliers on those rows are
  gone. An auxiliary that an equation outside the block shares stays a
  constraint. See :doc:`/ModelShapes/DSGE/Optimal policy`.

- **Changed: one smoother for the regime-switching filters, in which
  identities written in the model hold.** ``imm``, ``imm_o``, ``gpbn``
  (every history length), ``immn`` and ``crskfc`` now smooth pairs of
  linked regime histories against each pair's own prediction, with
  Kim's posterior weights. An identity the model states in every regime
  (a level and its growth rate) now holds in the smoothed path to
  rounding; before, it was off by up to about 1e-2. Likelihoods and
  filtered series are unchanged; all smoothed output of switching models
  changes. ``crskfc``'s default smoother is now this one. See
  :doc:`/WorkingWithAModel/Filtering`.

- **Fixed: on the planner-discount knife-edge, the commitment solution of
  a regime-switching model is the limit too, rest points included.** The
  continuation to the limit of the planner's discount started only when
  the solver reported unit-circle eigenvalues, which the switching
  solvers never do. It now finds the unit roots in each regime's solved
  transition of the state variables. Under ``{'m','scl'}`` and ``'bm'``
  the expansion points are extrapolated with the decision rules. On a
  switching tracker model the reporting variable's multiplier is now
  zero, and the solution equals the model without it.
  Constant-parameter models, and switching models off the knife-edge, are
  unchanged. See :doc:`/ModelShapes/DSGE/Optimal policy`.

- **New: ``@outside_optimization(home){B_F}``, a block that exempts
  only the players it names.** The parentheses, refused until now, name
  the players that are not constrained by the block's equations and
  take the variables in the braces as given; the other players keep the
  equations as constraints. The exempt player has no multiplier on those
  rows and no first-order condition with respect to the variables it
  takes as given: the game is the one written with the equations in
  ``@model`` and the variables left out of the exempt player's
  objective, under commitment, discretion and loose commitment, in Nash
  and leader-follower games and in ``@state`` games. A block naming
  every player is the every-player form ``@outside_optimization{W}``.
  The linear-quadratic solvers still refuse both forms. See
  :doc:`/ModelShapes/DSGE/Optimal policy`.

- **Changed: a constrained update that cycles between two branches takes
  the one consistent branch, and otherwise fails at once.** Where a
  constraint binds, ``imm`` iterates its update through the constrained
  step. When the update kept alternating between two branches of the
  constraints, it failed after 20 iterations (retcode 314). It now checks
  each branch's update against the constraints. When exactly one is
  consistent, the update takes it, and the period is listed in the
  filter's output field ``resolved_cycle_updates``. When neither or both
  are, the filter fails with retcode 314 at once. Wherever the update
  settled before, likelihoods and retcodes are unchanged bit for bit.

- **Faster: filtering under constraints.** Three overheads of the
  constraint engine are gone: the ``max()``/``min()`` of a model's
  constraints no longer go through name-value parsing, the projection
  onto the constraints no longer throws and catches an error to look up
  its tolerance, and a projection whose Gauss-Newton loop only returns to
  the same points stops at once instead of at its step cap. Likelihoods
  are unchanged bit for bit. On a switching NK-HANK with a ZLB and
  borrowing limits, an evaluation with the constraints enforced takes
  about 45% less time.

- **Fixed: ``utils.math_comp.max(A)`` and ``utils.math_comp.min(A)``
  with one argument** errored. They now reduce as
  ``max(A,[],'ComparisonMethod','real')``.

- **Fixed: discretion under automatic differentiation.** With
  ``'solve_derivatives_type','automatic'``, next period's discretionary
  policy entered the first-order conditions as a constant: its
  derivatives were lost, and order-2 solutions differed from the
  symbolic ones. Automatic and symbolic derivatives now give the same
  solution.

- **Fixed: models built at ``max_deriv_order`` 0.** Building a model
  without symbolic derivatives, to solve it with automatic derivatives,
  failed with an index out of bounds. It now builds and solves.

- **Changed: a plan and ``simul_regime`` that disagree about a regime
  are an error.** ``simul_regime`` used to overwrite, silently, the
  regimes a plan fixes (``append(plan, 'regime', t, r)``). Where the two
  give different regimes for the same period, ``simulate`` now stops
  with ``RISE:simulate:regimeConflict``, naming the period; where only
  one gives a regime, or they agree, it is used. The jump-off date is
  not checked. See :doc:`/WorkingWithAModel/Forecasting and simulation`.

- **Changed: a constrained update that fails outright fails the
  filter (retcode 314).** Where a constraint binds, the filter updates
  through the constrained step; when that update failed it fell back
  silently to the Kalman update, so the likelihood was not the
  constrained model's. It now stops with retcode 314
  (``CONSTRAINED_UPDATE_FAILED``), which an estimation treats as an
  inadmissible draw. An observed variable held at its bound and an
  ill-conditioned branch keep the Kalman update, as before. On the test
  models no update fails outright, and likelihoods are unchanged. See
  :doc:`/WorkingWithAModel/Filtering`.
- **Fixed: in ``@transition_functions`` a log variable means its level,
  as in ``@model``.** The transition matrix read a log variable as its
  log while the reported series (``zlb_tp_1_2``) read its level, so
  the two disagreed, and a transition function with a negative power of
  a log variable divided by zero at the steady state. Every block now
  reads a log variable as its level. **A model written for the old
  reading must be rewritten in levels** (``exp(theta*B)`` that stood
  for ``B^theta`` becomes ``B^theta``). See
  :doc:`/WorkingWithAModel/Time-varying transition probabilities`.

- **Fixed: under constraints, the regime-switching filters predict each
  regime with that regime's solution.** For each regime of the next
  period, the constrained prediction was taken in the most likely
  successor of that regime instead: the same regime when the transition
  probabilities are constant and persistent, often another one when they
  vary with the state or a regime stays with probability below 0.5. The
  prediction then disagreed with the covariance and the update of its
  branch, and the update through the constrained step fell back to the
  Kalman update (35 of 78 periods on a test HANK with time-varying
  probabilities; none now, and ten times faster). Constrained
  likelihoods of such models can change; models with constant,
  persistent probabilities are unchanged, bit for bit. Applies to every
  regime-switching filter, x-at-risk and the reduced models' filter. See
  :doc:`/WorkingWithAModel/Filtering`.

- **Changed: Taylor projection takes the expectation over next period's
  shocks.** The ``rtp`` paradigm set next period's shock to zero, so its
  solution had no risk correction at any order. The projection equations
  now hold in expectation (``rtp.shock_nodes``, a rule exact at the
  order): at the steady-state anchor order-2 Taylor projection is order-2
  perturbation, and re-anchored it carries the same precautionary drift.
  Under regime switching the regimes' policies now differ by their
  volatilities. ``'expectation','none'`` restores the old,
  certainty-equivalent projection; order 1 is unchanged. See
  :doc:`/WorkingWithAModel/Extending RISE through paradigms`.

- **Changed: under the OccBin paradigm a constrained variable lands on its
  bound when the binding branch is nonlinear in it.** The constant of a
  binding regime is taken where its constrained variables sit on their
  bounds, so a rule ``log(1+R) = gam*log(...)`` gives ``R = 0`` at the zero
  lower bound rather than the tangent value at the steady state. The new
  option ``'anchor_at_bound',false`` restores the classical OccBin constant.
  See :doc:`/ModelShapes/DSGE/Occasionally-binding constraints` and the
  tutorial ``ModelShapes/occbin_bound``.

- **New: equations outside the optimization.** An
  ``@outside_optimization{W}`` block holds equations that are solved with
  the model but are not constraints of the optimal-policy players, such as
  a welfare tracker kept for reporting: they carry no multiplier, and the
  variables in the braces, one per equation, are taken as given by the
  players. The policy problem is then the one without those equations.
  A variable taken as given that enters what a player optimizes is
  reported (``RISE:outside_optimization:feedback``), or refused with
  ``'outside_optimization_feedback','error'``. Exempting only some players
  and the linear-quadratic solvers are not available yet. See
  :doc:`/ModelShapes/DSGE/Optimal policy` (*Equations outside the
  optimization*) and the tutorial ``OptimalPolicy/OutsideOptimization``.
- **Changed: the multiplier charts take the exact trend exponents of the
  Lagrange multipliers; one numeraire per independent trend.** Under
  ``geometric_multipliers`` with ``multiplier_charts``, the exponents of the
  multipliers' trends come from the symmetry of the static system rather
  than from the solved growth factors, so exponents such as -1 or 0 are
  exact and zero or negative multipliers get theirs from their first-order
  conditions (``model_data.diagnostics.multiplier_trends``).
  ``multiplier_numeraire`` accepts a cellstr, one trending variable per
  independent trend, such as technology and the price level. The cases the
  numeraires do not span, multipliers that do not grow as a power law, and
  models where the exact route is declined are reported with warnings
  (``rise:multiplier_charts:notSpanned``, ``nonGeometric``,
  ``exactDeclined``). See :doc:`/ModelShapes/DSGE/Optimal policy` and the
  tutorial ``OptimalPolicy/NonstationaryRamsey`` (``howto_multiple_trends``).
- **Fixed: under discretion the mixed derivatives of next period's policy
  functions counted several times.** The first-order conditions of
  discretion and of stochastic replanning carry the derivatives of next
  period's policy; a mixed second derivative entered them twice (a mixed
  third one three or six times). Solutions at order 2 and above change
  wherever the policy functions take two inputs or more; first-order
  solutions do not. On the quasi-hyperbolic consumer of Krusell and Smith
  (2003), whose Markov-perfect equilibrium has a closed form, an order-2
  solution with a mixed derivative in the policy functions now matches the
  closed form to 5e-8, where its first-order coefficients were off by up to
  0.5. ``sstate_discr_policy_init`` now sets the initial level and first
  derivatives of the policy functions (it had no effect).
  See :doc:`/ModelShapes/DSGE/Optimal policy`.
- **Changed: every first-order solve goes through the model's configured
  solver; the linear-quadratic optimal-policy solvers state their own
  limits.** The self-consistent linearization solved the first orders of its
  fixed point with the general solver and refused ``loose_commitment`` and
  ``stochastic_replanning`` by name. It now uses the solver the model is
  configured with, whichever it is, as ``solve`` does
  (``rise.engine.dsge_tools.solve.configured_solver``). Each solver applies
  its own requirements: the linear-quadratic solvers need the lead
  derivatives of a pair of regimes, divided by the transition probability,
  to depend on the next regime only, and stop with
  ``RISE:looseCommitment:leadBlockDependsOnPreviousRegime`` (or
  ``RISE:stochasticReplanning:...``) when the matrices they are given do
  not have that structure; above order 1 they stop with
  ``RISE:looseCommitment:firstOrderOnly`` (or
  ``RISE:stochasticReplanning:firstOrderOnly``). On a linear-quadratic
  planner where the structure holds, loose commitment under ``{'m','scl'}``
  gives the general solver's answer under ``{'m','scl'}``. Nothing changes
  under ``'m'``.
  See :doc:`/ModelShapes/DSGE/Optimal policy`.
- **Changed: the self-consistent linearization and the Barthélemy-Marx
  anchors solve optimal-policy models.** ``{'m','scl'}`` and ``'bm'`` were
  refused on any model with an ``@optimization_problem``. A perturbation
  scheme does not depend on where the equations come from: the first-order
  conditions are now re-anchored like any other equations, under
  commitment, under discretion and wherever the commitment parameter
  switches, at every order, and with a common steady state both types
  return plain ``'m'``. Under discretion the policy functions of next
  period's regime, which the first-order conditions carry, are evaluated
  where the self-consistent linearization puts the leads, at the transition
  points; under ``'m'`` nothing changes.
  See :doc:`/ModelShapes/DSGE/Perturbation types` and
  :doc:`/ModelShapes/DSGE/Optimal policy`.
- **Fixed: Taylor projection's analytical Jacobian follows the lead point;
  a re-anchored root far from the policy is rejected.** The model's
  derivatives are evaluated where the next period's policy lands, which every
  coefficient moves; ``rtp.jacobian_analytical`` held them fixed (2.3% off the
  numerical Jacobian, and order-3 Newton stalled where a root exists). With
  the exact Jacobian, Newton beyond a fold of the order-2 equations in levels
  converged to absurd roots, so a converged re-solve is now kept only if
  today's policy at the new state stays within ``reanchor_max_jump`` (new
  option of ``rtp.paradigm``, default 0.5) of the last converged policy there;
  otherwise the step warns and reuses the last converged policy. In levels,
  order 3 now converges at every step from 0.9 to 0.3 of the steady-state
  capital. Written in logs (``@log_vars``), the growth model keeps the root of
  its projection equations everywhere, and order-2 Taylor projection is 40 to
  450 times more accurate than order-2 perturbation in logs. Write models in
  logs for Taylor projection far from the steady state.
  See :doc:`/WorkingWithAModel/Extending RISE through paradigms`.
- **Fixed: Taylor projection differentiates the next period's policy where
  the state lands.** The projection residual used the next period's policy
  and its derivatives at that policy's own anchor, not where the state lands
  one period later; the two coincide only at a steady-state anchor of a
  constant-parameter model. Re-anchored order-2 solutions were therefore less
  accurate than order-2 perturbation near the steady state (on a growth
  model started at 90% and 70% of its steady-state capital, largest
  consumption errors of 0.060% and 0.536% of the steady state against 0.006%
  and 0.188%); they are now more accurate (0.0033% and 0.138%). Order 3 is
  as accurate as order-3 perturbation or better from 90% down to 30%. From
  half the steady-state capital or less, the order-2 projection equations
  have no root near the policy: those re-solves warn and reuse the last
  converged policy, and order 3 is the order to use.
  See :doc:`/WorkingWithAModel/Extending RISE through paradigms`.
- **Changed: the RLL paradigm solves each local system with QZ and keeps its
  own constant by default.** The defaults of ``rll.paradigm`` are now
  ``'solver', 'rise_qz'`` (was ``'fi'``) and ``'use_steady_state', false``
  (was ``true``); both old values remain options. Functional iteration,
  capped at 100 iterations, failed at every re-linearization of the growth
  model far from the steady state, so RLL silently equaled first order
  apart from one warning; with ``use_steady_state`` true the step applied the
  local slopes around the global steady state, which was less accurate than
  first order. From capital at 0.9, 0.7, 0.5, 0.4 and 0.3 times its steady
  state, the largest error in consumption against perfect foresight (% of
  the steady state) is now 0.07, 0.73, 2.25, 3.28 and 4.27, against 0.11,
  1.16, 3.76, 5.94 and 9.05 at first order. At an occasionally binding bound
  the local system has no unique stable solution, so RLL keeps the previous
  local solution there, as before.
  See :doc:`/WorkingWithAModel/Extending RISE through paradigms`.
- **Fixed: each simulation under the RLL paradigm starts from the solved
  local policy.** ``rll.paradigm`` updated its stored solution in place, so a
  simulation started from the local policy the previous one ended with and
  the same simulation run twice could give two paths. Each simulation now
  works on its own copy, as under ``rtp.paradigm``.
  See :doc:`/WorkingWithAModel/Extending RISE through paradigms`.
- **Changed: a lead that a lag multiplies in a constraint of an optimal
  policy problem gets its own expectation variable.** In
  ``b{t} = b{t-1}/betta*(1+0.1*pie{t+1}) + ...`` the first-order condition
  for ``b`` would carry ``pie{t+2}``. RISE now gives ``pie{t+1}`` the
  variable ``LEAD_1_pie`` before it derives the conditions, as if
  ``EPIE{t} = pie{t+1}`` had been declared in the model. Solutions under
  commitment are unchanged up to terms of higher order; under discretion
  and stochastic replanning the solution is now that of the expectation
  form, which solves much faster. Terms in which a lag and a lead interact
  nonlinearly are left as written. See :doc:`/ModelShapes/DSGE/Optimal policy`
  ("A lead that a lag multiplies").
- **Fixed: re-anchored Taylor projection from a state far from the steady
  state.** Under ``rtp.paradigm`` (re-anchoring is the default), a growth
  model started at half its steady-state capital took a path on which
  consumption jumped up and capital ran down: the first re-solve started
  from the steady-state policy and converged to another root of the
  projection equations. Each re-solve now starts from the current policy
  evaluated at the new state, and only the states of the date-0 vector
  matter. The solved model also kept the anchors of the last simulation, so
  running the same simulation again could give another path; every
  simulation now starts from the solved polynomial.
  See :doc:`/WorkingWithAModel/Extending RISE through paradigms`.
- **Fixed: the heterogeneous-agent steady-state engine no longer runs
  away.** On a Huggett household with a borrowing limit of -2, the
  Anderson-accelerated time iteration of
  ``rise.engine.dsge_tools.sstate.het.engine`` could throw the policy to
  absurd values (steps of 4e35) and then sweep slowly without end. A mixed
  step much longer than the plain one is no longer taken, a sweep whose step
  blows up restarts from the best iterate, and repeated blow-ups stop the
  solve with the warning ``het:engine:diverged``.
  See :doc:`/ModelShapes/DSGE/Heterogeneous agents`.
- **Fixed: a plan with free shocks runs over its whole span on a
  regime-switching model.** On a regime-switching model, a plan whose shocks
  are free (``NaN``) holds nothing after the jump-off, the regimes being left
  to the chain, and its database lost those dates: ``simulate`` returned the
  jump-off alone, and a conditioned plan stopped after its last condition. The
  plan's database now keeps every date. Single-regime models were not
  affected. See :doc:`/WorkingWithAModel/Simulation plans`.
- **Fixed: simulate honors ``simul_history_end_date`` with a plain
  database, and forecast draws its shocks under
  ``forecast_shock_uncertainty``.** ``simulate`` ignored the date and
  started at the database's first date, pinned to the data after it.
  When the history runs past the database's first date, the simulation
  now jumps off at the date from the plan ``simplan.from_history``
  generates out of the history, as ``forecast`` does. A database that
  starts at the date is a conditioning database, as before. The shocks
  of such a plan are now drawn under shock uncertainty (the
  constraint-enforcing shocks excepted); before, they were pinned, so a
  forecast from a database gave the same path whatever
  ``forecast_shock_uncertainty`` said. ``simplan.from_history`` takes
  the regime and the shocks as optional arguments and no longer fails
  when the deepest lag predates the data.
  See :doc:`/WorkingWithAModel/Forecasting and simulation` ("Starting
  from history").
- **Fixed: simulate takes a recursive plan under a solve paradigm.** A plan
  built with ``'recursive', true`` carries a second page, and under a solve
  paradigm (for example the piecewise-linear paradigm on a model with a
  lower bound) ``simulate`` failed with a size mismatch in the sub-step
  forecast. A paradigm steps with the shocks of the model's own anticipation
  horizon: later pages that announce nothing are now dropped, so the
  recursive plan runs as surprises, and a page that announces something
  beyond the horizon is refused with
  ``RISE:simulation:planExceedsParadigmHorizon``. Under ``simulate``,
  page-1 entries are surprises whether the plan is recursive or not.
  See :doc:`/WorkingWithAModel/Surprises and announcements` ("Under
  simulate").
- **Fixed: the historical decomposition of a VAR gives the constant its
  own column.** ``historical_decomposition`` read the constant of a VAR
  from a place where the solver leaves zero, so the ``const`` column was
  zero and the structural-shock columns absorbed the constant (the
  columns still added up to the data). The columns are now ``init`` (the
  initial lags propagated), ``const``, one per deterministic variable,
  and one per structural shock, which contain the shocks' contributions
  only; they add up to the data. With regime switching the regimes are
  weighted by the smoothed regime probabilities and the output has one
  page instead of one per regime. ``residuals`` and ``structural_shocks``
  no longer contain the constant either. Results change for every VAR
  with a constant and for every switching VAR.
  See :doc:`/ModelShapes/Main Reduced form VAR Modeling`.
- **Fixed: tuning the independence sampler maximizes its acceptance rate.**
  With ``do_tuning``, ``imh`` used the random walk's rule (shrink the scale
  until the acceptance rate reaches 0.234), which collapses an independence
  proposal: on a multimodal target the scale fell to about 1e-200. It now
  searches the scale ``c0*2^k`` that maximizes the acceptance rate and fixes
  it; on Gaussian targets it lands within a fifth of the best scale. Warnings
  ``rsamplers:imhTuningFailed`` and ``rsamplers:imhTuningAtSearchEdge`` flag a
  search that cannot succeed. The random walk is unchanged.
  See :doc:`/Estimation/Posterior simulation` ("Tuning the proposal scale").

- **Changed: the Choleski identification of a reduced-form VAR follows
  the order in which the variables were declared.**
  ``identify(m, 'choleski', shock_names)`` and ``identification`` used to
  run the recursion in RISE's alphabetical order of the variables; they
  now follow the list given to ``rfvar_model`` (units, then variables, for
  ``prfvar_model``), as a recursive ``svar_model`` does. Shock *k*
  belongs to the *k*-th declared variable and takes the *k*-th name; the
  new option ``'ordering'`` sets another order. An unidentified
  reduced-form VAR uses the same recursive factor. Results change unless
  the variables were declared alphabetically; ``'ordering',
  sort(variables)`` reproduces the old recursion. Also, ``irf`` now labels
  the responses by shock slot, so it always agrees with the variance and
  historical decompositions, including with fewer names than variables.
  See :doc:`/ModelShapes/Main Reduced form VAR Modeling`.


- **Changed: template differentiation is the default.** The construction
  option ``template_diff`` is now ``true``: every model is differentiated one
  distinct expression at a time, and an equation with no repeated term is
  differentiated as before. Same derivatives (to rounding) at every order,
  with and without regime switching. ``dsge_model(file, 'template_diff',
  false)`` keeps the standard path.
  See :doc:`/ModelShapes/DSGE/Template differentiation`.

Release 20260831
----------------

- **New: template differentiation, construction option ``template_diff``.**
  Models that repeat a few equations -- the cells of a heterogeneous-agent
  model, the countries or sectors of a loop-generated model, or copies typed
  out by hand -- are differentiated one distinct expression at a time: RISE
  finds the terms that are the same expression up to a renaming of variables
  and parameters, differentiates each once, and evaluates its derivatives for
  all copies in one vectorized call. Nothing to annotate:
  ``dsge_model(file, 'template_diff', true)``. Same derivatives as before (to
  rounding) at every order, with regime switching (constant or time-varying
  probabilities) and in disk mode. On a two-asset HANK model with 150 cells,
  differentiation falls from 655 s to 40 s and peak memory from 29.8 GB to
  5.2 GB; Jacobian evaluation is 30 to 230 times faster across the models
  measured. Off by default.
  See :doc:`/ModelShapes/DSGE/Template differentiation`.
  Tests: ``models/dsge/template_diff``; tutorial:
  ``ModelShapes/template_differentiation``.

- **Fixed: options given inside the model cell no longer discard the
  name-value options.** ``dsge_model({file, name, value}, name2, value2)``
  rebuilt the construction options from the cell alone, silently dropping
  ``name2`` (e.g. ``'smooth_kinks', false``). The cell's names and values are
  now validated and merged with the name-value options; a runtime option in
  the cell, an unknown name, or the same option with two different values is
  an error. Options are meant to be passed as name-value arguments; cells
  holding the model line by line, several snippets or several files are
  unchanged.

- **New: ``rise.microfound`` derives Kuhn-Tucker (inequality) constraints.**
  ``.subject_to(irrev="I >= 0")`` used to be accepted and then silently
  dropped from the derived model -- you got a clean solution to the
  UNCONSTRAINED problem, with no error and no warning. Inequalities are
  now normalized to their ``g <= 0`` side, given a multiplier
  ``<agent>_lambda_<label>`` that enters the FOC of every control the
  constraint touches, and emitted as the complementarity row. Since the
  Lagrangian is ``L = U + lambda*g``, the derived multiplier is the
  negative of the textbook one, so ``lambda <= 0``, ``g <= 0``,
  ``lambda*g = 0`` is exactly ``max(lambda, g) = 0`` -- written out as that
  ``max``, which the parser already recognizes as a kink and routes
  through the occasionally-binding-constraint machinery. Both ``>=`` and
  ``<=`` are accepted (strict forms treated as non-strict). A slack
  constraint is inert: ``examples/microfound/rbc_irreversible.m``
  reproduces the plain RBC steady state with ``I = delta*k`` and a zero
  multiplier, while a binding floor pins the variable and returns a
  negative multiplier. Regression test
  ``models/dsge/microfound/testMicrofoundKKT``.

- **``rise.whatsnew(n)`` no longer spends a slot on an empty section.**
  Cutting a release leaves a freshly opened, empty ``Unreleased`` heading
  behind, and ``rise.whatsnew(1)`` used to print that -- the preamble and
  nothing else. Sections with no entries are now skipped, so ``n`` counts
  sections that actually carry entries.

- **Fixed/New: mixed OBC enforcement — a ``max``/``min`` kink and a
  regime-switching constraint can now coexist in one model, each enforced
  by its own mechanism.** Previously a kink on a switching model silently
  disabled the regime machinery: the constrained simulation path took the
  deterministic modal transition every period (a persistent chain never
  left its initial regime), a user ``simul_regime`` switch rule errored out
  with "Both simul_constraints and simul_regime were provided" (the kink
  had auto-set ``simul_constraints``), and user-declared ``simul_constraints``
  overwrote the auto-generated kink rows, leaving the kink enforced by
  nothing. The enforcement contract is now: **the regime comes first, the
  fudge shocks second**. Each simulated period, the regime-switching
  machinery selects the regime (stochastic transition draw per
  ``simul_stochastic_regimes``, or a switch rule reweighting the transition
  row), and the one-step forecast then resolves the *fudge-owned*
  restrictions within the selected regime — including inside the
  candidate-regime forecasts a switch rule compares, and under
  pre-specified regime paths. Ownership is per restriction row:
  kink-derived rows (with their auto ``FUDGE_FACTOR_*`` enforcers) are
  always fudge-owned; user rows are regime-owned when given in 2-column
  form on a switching model and fudge-owned otherwise; user
  ``simul_constraints`` on a kinky model now **merge** with the kink rows;
  the auto-compiled switch rule covers only the regime-owned rows; a user
  switch rule is accepted alongside fudge-owned rows and conflicts only
  with regime-owned rows (clear error). Pure kink-only and pure
  regime-switching-OBC models are unchanged (all locked regression
  numbers reproduce). Estimation/filtering violation-check forecasts
  still force the deterministic modal regime path for reproducibility.
  Worked example and regression test:
  ``models/dsge/obc/worked_examples/mixed_enforcement`` in the test suite.

- **New: ``threshold_chains``, a one-line declaration for threshold and
  smooth-transition VARs.** Available on every VAR constructor
  (``rfvar_model``, ``svar_model``, ``prfvar_model``, ``proxy_svar_model``):

  .. code-block:: matlab

     m = rfvar_model(endo, 'lag_length', 1, 'threshold_chains', struct( ...
             'name','sc', 'controlled_parameters',{{'b','c','covar'}}, ...
             'threshold_variable','ACR'));


  This was always expressible by hand — a two-state chain whose
  ``endogenous_probabilities`` are a steep logistic — and the machinery is
  unchanged: the sugar desugars into exactly that ``markov_chains`` entry
  before anything downstream sees it, so a sugared model reports the chain
  it actually has. What the sugar buys is that four things which are
  invisible when wrong are now impossible to get wrong: that the third
  argument of ``logistic`` is the threshold; that the two probabilities must
  be **complements** (a threshold rule is memoryless in ``s_{t-1}``, so two
  independent logistics silently give a smooth-transition chain, which is
  a different model); that ``delay = 1`` is written with **no lag operator**
  because the transition matrix built at ``t`` governs the regime at ``t+1``;
  and that a slope too flat for the scale of the variable leaves
  observations in the logistic's interior, where they become genuine
  regime *mixtures*.

  ``threshold_value`` and, for the smooth transitions, ``slope`` each take
  either ``'estimate'`` (the default — creates ``<chain>_bar`` /
  ``<chain>_gam``, free to be given a prior) or a **number**, which is
  written into the transition expression as a literal and creates no
  parameter at all. The asymmetry is deliberate: a fixed value that
  arrived as a parameter would be one more name the user must remember to
  ``set``, and forgetting is silent — the model builds, filters and returns
  a likelihood computed at whatever the uncalibrated value happened to be.
  That is the class of failure this sugar exists to remove, so a fixed
  threshold does not get a name.

  ``transition`` selects ``'hard'`` (default), ``'logistic'`` or
  ``'exponential'``; the hard case uses the new
  ``utils.smooth_transition.hard_threshold``, the ``g -> inf`` limit of the
  logistic floored into ``[1e-10, 1-1e-10]``, so there is no slope to
  choose. It is named ``hard_threshold`` and not ``threshold`` because
  MATLAB's Econometrics Toolbox ships a ``@threshold`` class and a class
  folder takes precedence over a path function — the shorter name is
  silently shadowed.

  Two implementation notes for anyone touching it. ``hard_threshold``
  deliberately does **not** wrap its comparison in ``double()``: RISE
  differentiates transition functions through ``@rsymbdiff``, which
  overloads ``gt`` but has no ``double``, so the cast would make the function
  unusable inside a ``@transition_functions`` block. And ``threshold_chains``
  is dropped from the options handed to the code generators
  (``rise.engine.parser.generator_config_pairs``), which declare exactly the
  configuration they accept.

  Covered by ``models/var/reduced_form/testThresholdChains.m`` (11 tests,
  including that the sugared chain reproduces the *lagged* rule, that the
  contemporaneous rule demonstrably does not fit, and that a model built
  with a fixed threshold gives the same likelihood as the same chain
  written out by hand). Documented in the reduced-form VAR chapter.

- **New: local projections, with state interaction and sign-restricted
  identified sets.** ``rise.engine.var.lptools.estimate`` runs the
  projections and ``rise.engine.var.lptools.identified_set`` searches
  rotations over them.

  A local projection is neither a VAR nor a univariate model, and it does
  not get a model class. There is no dynamic system, no state space and no
  likelihood: it is a stack of OLS regressions of leads on current and
  lagged variables, and nothing about it is sampled. What makes it
  structurally usable is that, with ``y_t`` on the right-hand side, the
  coefficient block on ``y_t`` at horizon ``k`` **is** the reduced-form
  impulse response matrix at that horizon (Plagborg-Møller & Wolf). So the
  estimator's whole job is to hand back ``{B_k}`` and ``Sigma`` — exactly the
  pair a rotation-based scheme consumes — and two functions are the honest
  shape for that.

  Restrictions are parsed by the same routine the VAR family uses
  (``vartools.read_identification_restrictions``, extracted from
  ``build_identification_engine`` where it had been a local function), so
  ``'GDP{0:2}@MP'`` means the same thing on both sides. ``identified_set``
  adds an ``elasticity`` bound on the ratio of two impact responses and a
  ``normalize`` option that rescales each accepted draw. Lag-structure
  restrictions are rejected: a projection has no lag polynomial.

  Two properties of the specification are enforced and documented rather
  than left to be discovered: horizon 0 is the identity by construction
  (its regression contains its own left-hand side), so ``Sigma`` comes from
  the ``k = 1`` residuals; and with a single pooled ``Sigma`` a
  state-dependent projection can only reveal state-dependence living in
  the **propagation** coefficients, never in the impact covariance. If a
  mechanism lives there, a switching VAR is the right tool.

  **Both local projections, not one.** Passing ``shock`` — the name of a
  *measured* shock series — switches ``estimate`` from the rotation-based
  form above to the Jordà (2005) / Ramey–Zubairy (2018) form, where the
  coefficient on that series **is** the structural response and no
  identification step exists. ``est.mode`` records which you ran. The two
  answer different questions and are not interchangeable, so they are
  labelled rather than blended: shock mode returns ``beta`` and ``beta_se``
  and leaves ``B`` empty, and ``identified_set`` rejects such an estimate
  outright rather than rotating nothing.

  Shock mode **drops the contemporaneous block**, and that is the whole
  design. ``s_t`` causes ``y_t``, so ``y_t`` mediates the effect being measured;
  conditioning on it closes the channel and drives ``beta`` toward zero — a
  bad control, and a silent one. Only lags of ``y`` are controls. A
  consequence worth having: ``beta`` at ``k = 0`` is a real estimate, the
  impact response, where lpvar mode's ``B{r}{1}`` is the identity and
  carries no information.

  Standard errors are Newey–West, truncated at ``k + 1`` by default
  (``hac_lags`` to override). Overlapping horizons leave the horizon-``k``
  residual serially correlated of order roughly ``k``, so OLS standard
  errors are wrong and get worse with the horizon; shipping ``beta`` without
  a HAC correction would have shipped an unusable number. lpvar mode
  reports no standard errors — its uncertainty is an identified set.

  **And the validation direction, with its precondition.** Running a
  projection on data simulated from a model you control — and checking it
  returns that model's own responses — is how you learn whether a
  specification can recover a mechanism before trusting it on real data.
  It works: on a known VAR(1) with a triangular impact matrix,
  ``chol(Sigma,'lower')`` recovers ``A`` to ``0.002`` and ``B_k*A`` matches the
  model at every horizon to ``0.004``.

  It works *only when the shock is recoverable from the observables
  handed over*, and the failure is silent. For
  ``y_t = e_t + 2 e_{t-1}`` — non-fundamental, with autocovariances
  identical to the invertible ``y_t = u_t + 0.5 u_{t-1}`` — the projection
  returns ``0.502`` where the truth is ``2``, reporting decay where the effect
  grows, at ``T = 200,000``. More data never helps. This is not a weakness
  of local projections: LP and VARs estimate the same object, so LP
  inherits invertibility problems in full and is not a robustness fix for
  non-fundamentalness. Supplying a shock series recovers ``[1, 2]`` on
  exactly that data — the fix is information, not a different estimator.

  Covered by ``models/var/local_projections/testLocalProjections.m`` (19
  tests, including that a linear projection on VAR(1) data reproduces
  ``B^k``, that shock mode recovers ``(B^k A)(:,1)`` exactly — impact vector
  included — that the model→data→projection round trip closes, that the
  non-fundamental IRF is *not* recoverable from history while shock mode
  gets it, that a normalization is exact, and that an elasticity bound
  binds). The bad-control test uses a *noisy proxy* deliberately: with the
  exact structural shock, ``e_t = inv(A)*(y_t - B*y_{t-1})`` makes the
  bad-control design singular, so the comparison would measure a
  minimum-norm artefact rather than the bias it claims to. Documented in
  the new *Local Projections* chapter; tutorial at
  ``rise-modern-tutorials/ModelShapes/local_projections``.

- **New: a policymaker can declare that it distrusts its model.** Add
  ``@robustness = rho`` to an ``@optimization_problem`` and the policy becomes
  robust in the sense of Hansen and Sargent: it is chosen to perform
  acceptably not only under the model as written, but under the least
  favourable model close to it.

      @objective   = -0.5*(pi^2 + lambda_y*y^2 + lambda_i*i^2),
      @robustness  = rho_m,
      @instrument  = i

  That is the whole of it. The objective stays the ordinary loss: **you do
  not write the adversary**, its instruments, or any entropy penalty.

  ``rho = 1/theta`` is the robustness intensity and it belongs to the
  policymaker in the way risk aversion does -- it prices how far the
  policymaker entertains its model being wrong. It must name a *parameter*,
  never a number, so it can be calibrated, estimated, and attached to a
  Markov chain, which lets a policymaker fear misspecification more in one
  regime than in another. **``rho = 0`` is expressed by omitting
  ``@robustness``**: a model that does not ask for robustness does not pay for
  it, no distortion variables are created, and the answer is the ordinary one
  rather than an approximation to it.

  **No new solver was needed.** A robust policymaker behaves as if an
  adversary were choosing the worst misspecification it can afford, which is
  a two-player zero-sum game, which the non-cooperative-games engine already
  solves. What was missing was not an algorithm but the writing-out, and that
  is what the parser now does: one distortion instrument per shock, named
  ``DISTORT_<SHOCK>``, entering *inside* that shock's own coefficient and at
  every occurrence of it; a one-period lag, so the adversary commits before
  seeing the innovation it distorts; the adversary's objective as the
  policymaker's, negated, with the entropy term carried by both; and the
  adversary one ``@order`` behind, since ``min_u max_w`` is the ordering the
  problem is written in.

  Each of those five changes the answer if it is wrong, and two do so
  silently -- the entropy term placed in one objective instead of both, and
  the adversary allowed to see the shock it distorts, are each worth about
  1e-2 in the policy rule while passing every check that does not involve
  robustness. The fourth cannot be got right by hand even in principle: a
  ``#`` definition may not contain variables, so a hand-written game must
  duplicate its objective and the two copies drift.

  The distortions are ordinary endogenous variables and appear in the
  solution, in impulse responses and in simulations. That is deliberate: with
  more than one shock the worst case is a *direction*, and which
  misspecification a policy guards against is usually the point of the
  exercise. On a small two-shock New Keynesian model the policymaker puts 88%
  of the distortion on the cost-push shock and 12% on demand -- not because
  demand shocks are smaller, but because the instrument can offset them.

  Robustness composes with everything already there: both solution routes,
  every timing protocol, and regime switching, including a switching ``rho``.
  It is bounded, though: past a model-specific frontier
  ``rho < 1/lambda_max(C'*Pbar*C)`` the adversary can buy unbounded damage and
  the problem has no value. The frontier is an eigenvalue condition, so it
  moves with the model; under switching it couples regimes; and it moves as
  it is approached, since it depends on the value under the robust policy
  itself. It can be small -- near 3e-4 on an estimated New Keynesian model --
  so bracket it rather than guess.

  Validated against a hand-written statement of the same game to 5e-14 over
  two multipliers and both timings, and on Chen-Kirsanova-Leith -- three
  named shocks, two Markov chains, four regimes -- to 4e-16. Independently,
  the underlying recursion agrees with a coupled robust Markov-switching
  Riccati, an algorithm that appears nowhere in RISE, to 1e-14.

  Refused, with a reason: ``@robustness`` on more than one player of a game.
  Each would fear a *different* misspecification, so the problem stops being
  zero-sum and becomes a fixed point in rules -- a different object, not a
  bigger version of this one.

- **Fixed: a switching SVAR with endogenous transition probabilities could
  not be built at all.** ``generate_svar_dsge`` appended the
  transition-probability names to the ``@endogenous`` list. The modern
  parser rejects TVTPs declared as endogenous variables, so construction
  stopped with *"The following transition-probability names appear in
  @endogenous ... must be declared in @transition_functions"*. There was
  no way round it from the user side: the names are generated, not
  written by hand. ``proxy_svar_model`` delegates to the same generator and
  was blocked identically.

  The same generator also emitted the equations as legacy ``! <eqtn>;``
  lines inside ``@model``, which still parse but raise
  ``rise:tvp_bang_prefix_deprecated`` on every build.

  Both are now the modern form — the names stay off ``@endogenous`` and the
  equations go in a ``@transition_functions`` block — matching
  ``generate_rfvar_dsge``, which had been migrated while the SVAR generator
  was left behind. Covered by
  ``models/var/structural/testSwitchingSvarTransitionFunctions.m``, which
  checks construction, that the names are not endogenous, that no
  deprecation warning fires, that the solve gives a regime-specific
  ``A0^{-1} diag(sigma)``, and that proxy-SVAR inherits the fix.

- **Fixed: an informative VAR prior could not be used on a switching
  RFVAR.** ``estim_var_prior`` accepts a ``+var_priors`` template (Minnesota,
  NIW, IW, Sims-Zha, ...). Every one of those templates already exposes
  ``evaluate_at_reduced_form(B, Sigma_u)``, and ``svar_evaluator`` already
  evaluates it at each regime's reduced form and sums — the standard
  reading of "the same prior applied in every regime". But only
  ``svar_model`` and ``proxy_svar_model`` registered a
  ``routines.var_prior_evaluator``, so an RFVAR never reached that path.

  Both ends then took the older single-regime branch. ``merge_var_priors``
  assigned the template's parameters under BARE names, where a switching
  parameter needs one entry per state (``<pname>_<chain>_<state>``), so a
  switching RFVAR died during prior setup with "parameter c_acr is not
  controlled by const"; and had it survived, ``log_prior_density`` would
  have called ``template.log_prior(param)``, which is dimensioned for a
  single regime's ``vec(B)``. Non-switching RFVARs were unaffected, which
  is why this stayed out of sight.

  ``rfvar_model`` now registers
  ``+rise/+engine/+var/+var_priors/rfvar_evaluator.m``, the RFVAR
  counterpart of ``svar_evaluator`` and the simpler of the two: a
  reduced-form VAR is parameterized by its reduced form, so ``B`` is read
  straight off the ``c_*`` / ``b<lag>_*`` parameters with no ``A0`` to invert,
  and ``Sigma_u`` off ``covar_*`` (via the shared
  ``vartools.rfvar_sigma``). The evaluator branch of ``merge_var_priors``
  now also emits the RFVAR-shaped names (``b<lag>_*``, ``covar_*``) and
  keeps taking start values from the template's ``mean`` where it carries
  the parameter, so start values do not regress. Covered by
  ``models/var/reduced_form/testReducedFormVarPrior.m``, which checks the
  density against a hand-built two-regime sum.

- **Fixed: a reduced-form VAR's residual covariance never reached its state
  space, so every likelihood-based RFVAR route ran with ``Sigma = I``.** An
  RFVAR carries one residual per equation, correlated across equations,
  parameterised by ``covar_sigma_<i>`` (variances) and ``covar_<i>_<j>``
  (covariances). But ``generate_rfvar_dsge`` writes each residual into the
  generated model with a unit loading (``... + e1{t}``), so those parameters
  appear in no equation, and nothing read them back. The derivative solve
  therefore returned a shock block equal to the identity no matter what the
  parameters said.

  The consequences were silent and wide. The likelihood was *exactly* flat
  in ``covar_*`` — moving a variance from 0.04 to 16 left the log-likelihood
  unchanged to the last digit — so any ML or Bayesian estimation of those
  parameters was estimating on a flat surface and returning its seed.
  Regime-dependent covariances (``controlled_parameters`` containing ``'s'`` or
  ``'covar'``) were a no-op, which made regime heteroskedasticity in an RFVAR
  a no-op. Identification read ``Sigma_u = I`` and rotated the wrong object.
  Only the OLS short-circuit escaped, because ``estimate_var_ols`` bypasses
  the parameter dictionary and writes Sigma into the state space by hand.

  Fixed by ``+rise/+engine/+var/+vartools/apply_reduced_form_covariance.m``,
  called from ``dsge_tools.solve.solve`` after a successful solve and before
  identification: it assembles each regime's Sigma from the ``covar_*``
  parameters and writes a Cholesky factor into that regime's structural
  shock columns, the same block ``apply_identification`` later rotates. SVAR
  and proxy-SVAR are untouched — they already build their impact block
  structurally as ``A0^{-1} diag(sigma)``. Panel VARs are covered too:
  ``generate_prfvar_dsge`` expands the endogenous list to ``<unit>_<var>`` and
  delegates to the RFVAR generator, so their ``covar_*`` names resolve
  against ``var_block_names`` like any other RFVAR's.

  **No behaviour change for an uncalibrated covariance.** A ``covar_*``
  that was never given a value (NaN) keeps its identity entry — unit
  variance, zero covariance — which is precisely what the generated model
  hardwired before any of this was transmitted. Forcing people to
  calibrate a covariance matrix by hand to get a VAR to solve would be a
  tax on everyone to catch a mistake almost nobody makes; ``estimate``
  seeds ``covar_*`` from an OLS pre-fit anyway, so the only affected
  population is hand-built models, and those keep working exactly as
  before. NaN means "unspecified", not "invalid": a value that *was*
  supplied and cannot be a covariance (negative variance, an inconsistent
  pair, Inf) still stops the solve, with the new retcode
  ``RFVAR_COV_NOT_PD`` (264).

  Covered by ``models/var/reduced_form/testReducedFormCovariance.m``
  (per-regime Sigma, the panel case, identity defaulting for unset and
  partially-set covariances, rejection of a supplied-but-invalid one,
  likelihood responsiveness, and the shorthand forms).

- **Fixed: ``controlled_parameters`` entries naming a VAR covariance were
  silently ignored.** In ``assign_markov_parameters``, the shorthand expander
  tested ``startsWith(param,'c')`` before ``startsWith(param,'covar')``. Since
  ``'covar'`` also starts with ``'c'``, every covariance token fell into the
  constant/deterministic branch, matched neither ``'c'`` nor ``'c('``, and came
  back as an empty expansion — which the caller could not distinguish from
  a shorthand that legitimately matched nothing, so the parameter was left
  uncontrolled without a word. This made the entire ``'covar'`` family
  (``'covar'``, ``'covar(diag)'``, ``'covar(offdiag)'``, ``'covar(y1,y2)'``) dead
  code, and it also swallowed a fully typed-out name such as
  ``'covar_y1_y2'``, so an off-diagonal could not be made to switch by any
  route. ``'s'`` was unaffected — it is matched later in the chain — which is
  why the diagonal appeared to work.

  The covariance branch now precedes the constant branch, a typed-out
  ``covar_*`` name resolves to itself, and ``'covar(y1,y2)'`` strips its
  parentheses correctly (the old slice kept a leading ``(``). Related: a
  ``controlled_parameters`` entry that expands to nothing is now an error
  (``assign_markov_parameters:emptyExpansion``) instead of a silent drop —
  the check is on the expansion before the existence filter, so a shorthand
  may still generate candidates that do not exist, as the normalised ``a0``
  diagonal does. Covered by
  ``models/var/reduced_form/testReducedFormCovariance.m``.

- **Fixed: a switching parameter inside a player's objective was frozen at
  its first-regime value by the linear-quadratic game solver.** In a
  non-cooperative game solved with ``solver = 'loose_commitment'`` or
  ``solver = 'stochastic_replanning'``, each player's objective curvature was
  read once, at regime 1, and reused in every regime. The constraint blocks
  were always read at the regime they belong to; only the objective rows
  were not. A model whose loss weights switch -- a policymaker who penalises
  an instrument more in one regime than in another, say -- was therefore
  solved against the wrong loss in every regime but the first, and the
  general route (``solver = 'mfi'``) disagreed with it.

  The answer was wrong rather than merely imprecise: it was not an
  equilibrium of the game, and the wronged player's own best response
  departed from it. The error grows with how much the objective switches and
  vanishes when it does not, which is what kept the condition out of sight.
  Models with regime-invariant objectives are unaffected and no existing
  result changes.

  Fixed in ``+rise/+engine/+dsge_tools/+optimal_policy/lq_game_engine.m`` by
  carrying the regime on the objective block, as the constraint blocks
  already did. Covered by ``rf3sw_obj.rs``, a game whose loss weights switch
  and whose ``@model`` block does not.

- **New: the linear-quadratic optimal-policy solvers now accept
  non-cooperative games.** ``solver = 'loose_commitment'`` and
  ``solver = 'stochastic_replanning'`` previously errored with "cannot solve
  non-cooperative games"; both now handle any number of players at the same
  ``@order`` (same-level Nash), any number of Markov regimes, forward- or
  backward-looking constraints, and either equilibrium concept
  (``solve_policy_equilibrium`` ``'MPE'`` or ``'OLE'``). No new solver name was
  introduced: the re-optimisation protocol is the solver, the number of
  players is a property of the model, and the two are orthogonal.

  The point of the specialisation is that **it never evaluates a
  policy-derivative (``f1942``) function**. The general route approximates
  every player's decision rule by a polynomial and iterates on that
  approximation, re-evaluating the model's derivatives and re-solving the
  steady state at each pass. In a linear-quadratic problem those rules are
  linear and the coefficients being searched for are entries of the solution
  matrix, so the whole apparatus reduces to two blocks written in closed
  form: the rivals' reactions entering an augmented lag block, and the total
  derivative of the constraint that a re-optimising player prices. Those two
  are held fixed while an ordinary constant-coefficient
  rational-expectations problem is solved, then updated -- the same
  two-level structure the general route uses, and what selects the
  mean-square-stable root.

  On a two-country central-bank game with 34 endogenous variables and two
  regimes, under discretion with Markov-perfect play, the specialised route
  is 8-9 times faster (4.4-4.9 s against 39-40 s) and agrees to 1.5e-06,
  which is the *general* route's own accuracy. It also solves at default
  settings where the general route returns retcode 210 and needs
  ``solve_discretion_maxiter`` lifted from 100 to 500. On small models the
  advantage disappears: warm, the two are within a factor of two of each
  other at 13-17 equations, because fixed costs dominate.

  New engine at
  ``+rise/+engine/+dsge_tools/+optimal_policy/lq_game_engine.m``, with the
  model-ordered blocks in ``lq_game_blocks.m``. The single-planner path is
  untouched -- a one-player model never enters the game branch, and the
  existing 26 optimal-policy regression tests pass unchanged. The parser now
  records the per-player bookkeeping the engine needs (``player_wrt``,
  ``player_instruments``, ``player_order``) instead of leaving it to be
  recovered from equation-tag strings, and ``set_lq_game`` replaces the
  blanket refusal in ``set_loose_commitment``.

  **No new options.** The coefficient loop is the same loop, over the same
  objects, as the general route's, so it reads the same two settings:
  ``solve_discretion_maxiter`` for its iteration cap (floored at 50, exactly as
  in ``run_loop_local``) and ``solver{2}.TolFun`` for its tolerance. Tighten
  either and both routes respond identically. The one genuinely
  engine-specific choice -- which solver handles the inner
  constant-coefficient problem, since the name that selected this engine
  would not dispatch there -- travels in the solver's own option struct
  rather than becoming a property of the toolbox:
  ``solve(m,'solver',{'stochastic_replanning','InnerSolver','mnk'})``,
  defaulting to ``'mfi'``.

  One case is refused with a reason rather than answered. **Loose commitment
  with a commitment probability strictly inside (0,1)**, for a game: the
  probability would have to scale each player's forward multiplier block,
  which cannot be expressed in the form the shared first-order engine reads
  -- use ``stochastic_replanning``.

- **Settled: the promise coefficient is an EXPECTATION, not the realised
  value — and the last refusal is gone.** The linear-quadratic game engine
  used the promise block at the realised regime. The parser uses its
  expectation over the successor regime, conditional on landing in one that
  honours the promise:

      Abar(r_t) = sum_j Q(r_t,j) 1{gate(j)>0} A(j)
                  / sum_j Q(r_t,j) 1{gate(j)>0}

  The two readings coincide unless the promise block VARIES across the regimes
  that honour it, which needs three things at once: a hierarchy (so a leader
  holds promises on its followers' conditions), a commitment indicator that
  differs across regimes (so a promise can span a change of status), and a
  switching structural coefficient (so its value differs by regime). Any two
  agree; the triple did not, by 1.3e-03.

  The relationship was verified in closed form before anything was changed:
  with the structural chain at (0.9/0.1, 0.2/0.8) inside the commitment chain,
  the parser's coefficients are exactly that conditional mean of the engine's,
  to six figures in both commitment regimes. The engine now reads it the same
  way and the two routes agree to 7.3e-09. Fixture ``game_fllfsw``.

  Because the transformation is the identity when the block does not vary
  across honouring regimes, nothing that already agreed moved -- which is also
  why the case had to be constructed before the difference could be seen.

  NOTE for whoever revisits this: it makes the games engine AVERAGE where the
  single-planner promise block was fixed to use the REALISED block (see the
  ``sum_j Aplus{r0,j}`` defect). The objects differ -- this average is over the
  successor conditional on honouring, that one was over where the economy
  might land when r_t was already known -- but they rhyme, and the general
  route was taken as authoritative here rather than derived from first
  principles. One nested function (``promise_expectation``) to revert.

- **New: three or more ``@order`` levels are covered.** Every hierarchy fixture
  had exactly two. The construction is a forward pass over processing order and
  is generic in the count -- at the third level the row map resolves rows
  belonging to two different earlier groups rather than one. ``game_rf3x``
  (orders 1/2/3, one player each) prices 3 / 7 / 15 rows and agrees with the
  general route at 2.1e-10 (MPE) and 4.3e-10 (OLE). Checked rather than assumed,
  because "generic by construction" was also true of the Markov-perfect
  cross-chain before ``game_fl3lf`` disproved it.

- **Settled: the re-optimiser faces the RESET rule, and the last refusal is
  lifted.** Two channels in the linear-quadratic game engine stand in for what
  the general route evaluates with ``f1942``: the re-optimiser's continuation,
  and the Markov-perfect reaction of a rival's instrument. That machinery
  exists for DISCRETION and takes the non-multiplier variables as its inputs by
  construction -- a planner that re-optimises has discarded the promise it
  inherited, so the continuation it faces is a function of the genuine states
  alone.

  That is the reset projector ``R = diag(I,0)`` the engine already owns:
  ``loose_commitmentize`` applies exactly it, but only AFTER the fixed point.
  Inside the coefficient loop both channels were substituting the RAW rule.
  They now substitute ``H*R``.

  It is a no-op almost everywhere, which is why it went unseen. Under pure
  discretion ``H`` already equals ``H*R`` -- the multiplier columns of the rule are
  zero at the fixed point. Under pure commitment the discretion gate is zero
  and the continuation term is never used. Only a MIXED indicator reaches the
  case that matters: a discretion regime substituting forward into a COMMITTED
  successor, whose rule does respond to the promises made today (those columns
  measure 2.57, not zero).

  It is also read only by a player whose differentiation targets include a
  multiplier. Own-level multipliers -- a player's own and its same-``@order``
  peers' -- are never targets, so only a player above the bottom level reads
  them, and for a single-level Nash game the whole thing is a no-op.

  Costs, measured on a forward-looking two-level game with a regime-varying
  indicator: 0.37 of a leader's contemporaneous block from the continuation
  channel, and 1.5e-01 of the solution from the cross-chain. With both reset,
  the routes agree to 5.3e-09 (two players) and 2.2e-08 (three), and the
  refusal of a hierarchy with a regime-varying commitment indicator is lifted.

  New fixture ``game_fl3lf``: two leaders sharing a level above one follower,
  forward-looking. It is the only shape in which both channels are live on the
  same player -- a two-player hierarchy leaves the leader's rival set empty,
  and ``game_rf3lf`` has the two leaders but is backward-looking with zero
  multiplier columns. The cross-chain defect was invisible until it existed.

  Nothing in the general route changed, and nothing needed to.

- **Fixed (BOTH routes): the commitment indicator was mis-dated in a leader's
  forward block.** It is one symbol doing two jobs, and the right date depends
  on which row it sits in.

  On a STRUCTURAL row it asks "does the planner honour an inherited promise
  NOW", which is a property of ``t``, and ``date_switching_parameters`` has always
  -- correctly -- left it there. On a row that is a previously eliminated
  player's FIRST-ORDER CONDITION, which exists only under a hierarchy, it asks
  something else: that condition carries the FOLLOWER's promise term, and the
  block it sits in is multiplied by ``lambda{t+1}``, so what a leader prices is
  ``beta*E_t[lambda^L_{t+1} * c_{t+1}]`` -- "will my follower, NEXT period,
  honour the promise it is making now". The indicator belongs to the regime
  realised at ``t+1``.

  Left undated it was read at ``r_t``, and being a bare parameter it carried no
  tag for the forward shift to move: the same mechanism as the ``dc/du^rival``
  mis-dating. Worth 0.5 of the leader's lead block on a forward-looking
  two-level game, and it applied to the GENERAL route as much as to the
  linear-quadratic one -- any leader-follower model with a switching commitment
  indicator was affected whichever solver was chosen.

  The dating pass now splits by row. That is what makes the fix safe: a single
  planner and a same-level Nash game have no rows above the structural block,
  so neither can be reached at all, by construction rather than by testing.
  Measured on a forward-looking hierarchy under commitment: the assembled
  system and the parser's went from 9.1e-03 apart to 4.9e-15.

- **Hardened: the date shift refuses to run on an undated switching
  parameter.** Every coefficient that has to travel from ``t`` to ``t+1`` in the
  optimal-policy pipeline goes through one ``regexprep`` in the parser,
  ``local_lead``, and that expression moves *time-tagged* atoms only. A bare
  parameter carries no tag, so it does not move: it stays bound to the
  equation's own regime while every coefficient beside it lands at ``r_{t+1}``.
  For a constant parameter that is harmless. For a switching one it is a
  mis-dating, and an invisible one -- nothing errors, nothing warns, and it
  cannot be seen at all with a single regime. Two independent defects have now
  come through that funnel.

  ``local_lead`` now errors, naming the parameter, if it is handed a block still
  carrying an undated switching symbol. Exactly two are exempt, because their
  timing genuinely is not the constraint's: ``commitment``, which asks whether an
  inherited promise is honoured *now*, and ``ole``, which selects the equilibrium
  concept and is the same at every date. Models that were correct are
  unaffected -- the guard is silent once the dating pass has run.

- **New: leader-follower games on the linear-quadratic route.** Players with
  distinct ``@order`` values -- a Stackelberg hierarchy rather than Nash -- were
  refused by ``set_lq_game``, on the grounds that a leader differentiates its
  followers' conditions and so brings in derivatives of every tracked decision
  rule, and second derivatives of them. Under LQ those rules are linear and
  their coefficients constant, so the second derivatives vanish; and the
  followers' conditions are rows of the very system the engine writes, so a
  leader's dated blocks are a **row-restriction of the matrices being built**.
  Assemble the players in processing order -- ``@order`` descending, which is the
  order the parser eliminates them in, hence **lower ``@order`` = leader** -- and
  the rows a leader needs are already there. No new algebra, and it composes
  with the timing exactly as Nash does.

  Verified block by block against the parser at its converged coefficients on a
  two-player hierarchy: for the follower ``A0``/``A+``/``A-`` differ by
  ``0``/``6.3e-11``/``0``, for the leader by ``9.1e-12``/``0``/``0``. Against the general
  route the solutions agree to ``1.6e-09`` backward-looking, ``1.9e-10`` with a
  switching coefficient on a rival's instrument, ``1.4e-07`` forward-looking
  (the reference's own accuracy), and exactly under commitment.

  Beware the indexing: ``player_order``, ``player_constraint_rows``,
  ``player_multipliers`` and ``mult_groups`` are in *processing* order while the
  parser's declaration list is not. With one player per level the two coincide;
  ``game_rf3lf`` in the regression suite is the fixture where they do not.

  Refused, and deliberately not papered over: **a hierarchy whose commitment
  indicator differs across regimes**. A leader prices its followers'
  conditions, and a follower's condition carries a promise term gated by that
  indicator. Read as a constraint at ``t+1`` -- which is what the leader's
  forward block is -- the gate belongs to the regime realised at ``t+1``; the
  parser reads it at ``r_t``, because the indicator is a bare parameter with no
  time tag and ``date_switching_parameters`` excludes it *by name*, deliberately
  (in the single-planner promise block, "does the planner honour an inherited
  promise **now**" really is a question about ``t``). The readings coincide while
  the indicator is constant across regimes -- pure commitment and pure
  discretion solve normally -- and part company by ``0.5`` of the leader's lead
  block as soon as it is not. Making the engine agree by copying a dating it
  does not believe is precisely what hid the mis-dating below for as long as it
  hid it, so it refuses instead and says why.

- **Fixed: ``MPE`` no longer differs from ``OLE`` under commitment on the
  linear-quadratic game route.** Markov perfection is a discretionary concept
  and the general route says so by construction: it enters its coefficient loop
  on ``any(commitment ~= 1)`` and never on the equilibrium concept, so with a
  commitment regime everywhere the tracked policy derivatives stay at zero and
  the cross-player chain term is absent. The linear-quadratic engine computes
  those derivatives itself, so nothing stopped it applying the term anyway --
  and it did, putting it ``1.6e-02`` of the decision rule away from the general
  route on a two-player backward-looking game under ``solve_policy_type = 'ramsey'`` with ``'MPE'``. Both routes now agree exactly there.

- **Fixed: the linear-quadratic game engine wrote its first-order-condition
  rows negated relative to the parser.** A condition and its negation are the
  same equation, so for a Nash game the choice was free and invisible. It is
  not free once a leader prices its followers' conditions: the leader's row
  then mixed structural-multiplier columns carrying one sign with
  follower-multiplier columns carrying the other. Every condition is now
  written the way the parser writes it, which also means the two routes report
  the same shadow prices rather than a sign flip on the rows a leader carries.

- **Fixed (tests): two refusal tests in the linear-quadratic game suite never
  reached the solver.** They passed a handle of the form
  ``@()evalc('solve(m,...)')``; an anonymous function captures only the variables
  its body *mentions*, and a name inside a string literal is not mentioned, so
  the handle raised "Unrecognized variable m" wherever it was called and the
  "did it refuse?" assertion was satisfied by the wrong error. Both are now
  real.

- **Fixed (BOTH routes): the Markov-perfect chain term was mis-dated when a
  rival's instrument carries a switching coefficient.** Player ``p``'s
  optimality condition carries ``beta*E_t[lambda^p_{t+1}' * dc(t+1)/dy_t]``,
  and ``c(t+1)`` is the constraint block at regime ``r_{t+1}``. Its lag block and
  its derivative with respect to a rival's instrument are coefficients of the
  *same* date-``t+1`` equation, so both belong at ``r_{t+1}``. Only the first was.

  In the parser the cause was mechanical: ``mpe_envelope_inject``
  (``create_optimal_policy_first_order_conditions.m``) builds the cross term at
  ``t`` and shifts it forward with ``local_lead``, which rewrites *time-tagged*
  atoms only -- the multiplier and the variables. A bare parameter carries no
  tag, so the shift walked past it and it stayed bound to the equation's own
  regime. The fix tags it at ``t`` before the shift, so the shift carries it to
  ``t+1`` with everything else.

  Measured against the defining property of a Markov-perfect Nash
  equilibrium -- each player's rule must be its best response to the rivals,
  which with the rivals substituted into the law of motion and the loss is an
  ordinary single-player Markov-switching regulator -- on a two-player
  reduced-form switching game (exact answer 6.3e-12):

  .. list-table::
     :header-rows: 1

     * - route
       - before
       - after
     * - linear-quadratic
       - 4.3e-05
       - 5.5e-12
     * - general (``mfi``)
       - 4.3e-05
       - 1.9e-12


  and the two agree to 2.1e-10. When the rival's coefficient does not switch
  the two readings coincide and nothing changes.

  **Cost.** Only SWITCHING parameters are tagged, and only those appearing in
  ``dc/du^rival``, so the fix mints at most one auxiliary variable per such
  parameter and nothing otherwise. Measured: ``rf3`` (13 vars) unchanged,
  ``game_rf3p`` (18, three players) unchanged, ``game_2country`` (34, two
  regimes, switching Phillips slopes) unchanged; only the two fixtures whose
  rival coefficient actually switches grow, by 1 and 2 variables. Parse and
  solve times unchanged within noise.

  This also corrects an earlier entry which described the ``r_t`` reading as
  "settled by measurement". What had been measured was which dating
  reproduced the general route, on the assumption that the general route was
  the reference -- the assumption that needed testing. The equilibrium
  residual test is in the regression suite as ``mpe_br_residual.m``.

  Unrelated but worth recording, since it was conflated with the above during
  diagnosis: ``rise_1`` is not designed for regime switching, especially with
  switching in the dynamics or non-extreme transition probabilities, and
  scores 2.6e-03 on this test for that reason alone.

- **Fixed: the mean-square-stability check was dead on the grand replanning
  loop.** ``+rise/+engine/+dsge_tools/+solve/solve.m`` switches
  ``solve_check_stability`` off for the duration of that loop and restores it
  afterwards, but the post-loop call to ``is_stable_system`` sat *between* the
  loop and the restore -- and ``is_stable_system_engine`` returns "stable"
  immediately whenever that flag is false. The check was therefore inert on
  every discretion, loose-commitment and stochastic-replanning solve routed
  through the loop, which reported retcode 0 on solutions its own test
  rejects.

  Found by diffing the two routes on a two-player game with a commitment
  regime. Both returned the **same** solution -- 6.0e-08 apart, which is the
  general route's own accuracy, and each satisfies the system defined by the
  other's coefficients to 8e-09, so there was no multiplicity to choose
  between -- and ``is_stable_system`` rejected it with a mean-square-stability
  radius of 2.19. Only the linear-quadratic route, which is checked on the
  ordinary path, said so.

  The restore now happens before the call. **This is a behaviour change:**
  models previously reported as solved on that path may now return retcode

  25. That is the correct answer -- the model has no mean-square-stable
  equilibrium at that calibration -- but it will look like a regression to
  anyone who has been relying on the silent pass.

- **New: ``solver = 'stochastic_replanning'``, a linear-quadratic solver for
  the stochastic-replanning problem**
  (``+rise/+engine/+dsge_tools/+optimal_policy/replanning_solver_h.m``). It is
  the loose-commitment engine with the scalar ``gamma`` replaced by a
  per-regime commitment indicator read off the same ``commitment`` parameter,
  which changes three things. The transition matrix is unrestricted:
  ``loose_commitment`` requires ``Pr(commit next)`` to be identical from every
  regime and errors otherwise, because ``gamma`` is one number, so persistent
  and asymmetric replanning processes were previously out of reach. The
  commitment chain also need not be the only chain: everything is indexed by
  the composite regime, so any number of other Markov processes may switch
  the structural matrices, the loss weights or the shock loadings alongside
  it, with any number of states each. The reset lives in the equations
  rather than in a post-solve projection — ``GAM1(y,lamb)`` is gated by the
  indicator, so a replanning regime's zero multiplier response falls out of
  the fixed point and the reported ``Tz`` *is* the fixed point, unlike under
  ``loose_commitment``. And it remains available under ``no_u_turn``, since like
  the loose-commitment engine it computes the derivative of the future policy
  function inside its own fixed point and never needs the discretionary
  (``f1942``) machinery that ``@no_u_turn = true`` drops at parse time; on a
  model parsed that way with a persistent chain it is the only route that
  solves at all, ``loose_commitment`` refusing the chain and the general route
  refusing the missing derivatives.

  Validated against the general (parser + MSRE) route over the full reported
  solution, state columns and shock impacts alike: agreement to 1e-8 or
  better at ``gamma = 1``, at ``gamma = 0``, on i.i.d. chains with
  ``0 < gamma < 1`` — where ``loose_commitment`` disagrees with the general route
  by O(1), because it is solving the other problem — and on persistent chains
  ``loose_commitment`` cannot express. Indicatively around 50x faster than the
  general route on the test fixture, though that model is far too small to
  say anything about scaling.

  ``loose_commitment`` is untouched and reproduces its previous solutions
  exactly. Weighting the policy-derivative channel by the probability that
  the *successor* replans, rather than by the current regime indicator, turns
  the new solver back into loose commitment — exactly, on the i.i.d. edge.
  That variant is deliberately not exposed, because off that edge it has
  nothing to be validated against.

- **Fixed: the loose-commitment promise term used an average of forward
  matrices where the realised one is required**
  (``+rise/+engine/+dsge_tools/+optimal_policy/solver_h.m``). ``lambda_{t-1}``
  enters the date-``t`` optimality condition because ``y_t`` sat in the lead
  term of the date-``(t-1)`` constraint, so the promise term carries the
  forward block that constraint attached to ``y_t``, de-weighted:
  ``Aplus{iota,r_t}/q(iota,r_t)`` with ``iota = r_{t-1}``. The engine instead
  formed ``sum_j Aplus{r_t,j}``, marked in the source as an approximation.
  The two coincide whenever ``A+`` does not switch, which is why
  constant-coefficient loose commitment was never affected; they do not
  coincide once it switches. On a model with a lead-dated switching
  forward coefficient (``a{t+1}*x{t+1}``, i.e. ``A+_{ij} = q_ij*A+_j``) the
  decision rules were off by about six percent against the general
  solver, at both ``gamma = 1`` and ``gamma < 1``; they now agree to solver
  tolerance. The required block is now computed from the assembled
  derivatives and, where it cannot be represented, the engine says so:
  with a current-dated switching coefficient (``a*x{t+1}``, the default
  reading, giving ``A+_{ij} = q_ij*A+_i``) the promise matrix is indexed by
  the regime in which the promise was *made*, which a solution indexed by
  the current regime alone cannot carry; the historical approximation is
  kept and a warning names the general solver as the alternative. Models
  whose ``A+`` does not switch, including both existing loose-commitment
  fixtures, reproduce their previous solutions bit for bit.

- **Documented: loose commitment and stochastic replanning are different
  problems.** They are not two solvers for one problem, and at
  ``0 < gamma < 1`` they do not agree. Under loose commitment the
  re-optimization probability ``gamma`` is a scalar and re-optimization is
  not a regime but a state value (``lambda_{t-1} = 0``), so one decision
  rule per structural regime suffices; under stochastic replanning
  ``commitment`` is a regime indicator with its own rule, and the chain may
  be persistent, endogenous or time-varying. Both carry the same
  future-policy derivative — the ``f1942`` hyperparameters are filled from
  the solution's own ``Tz`` — and differ only in the scalar in front of it,
  ``1-gamma`` against ``1-commitment_t``. Those coincide only at ``gamma = 0``
  and ``gamma = 1``, where the two routes agree exactly and each solution
  satisfies the other's equilibrium conditions. See *Optimal policy* in
  the documentation for the full comparison and the measurements.

- **Fixed: an empty shock horizon no longer collapses ``siz.nz``**
  (``+rise/+engine/+dsge_tools/+solve/rehash_topology.m``). A model with no
  exogenous shocks has an empty ``exogenous.shock_horizon``, so
  ``max(shock_horizon(:))`` returns ``[]`` rather than ``0``. Since scalar plus
  empty is empty in MATLAB, the augmented-state size
  ``nz = np+nb+1+ne*(1+shock_horizon)`` silently evaluated to a 0x1 empty
  instead of a size, and so did anything computed from it downstream. The
  horizon is now normalized to a scalar zero before it enters any
  arithmetic, matching what ``load_stochastic_solution`` and
  ``set_simulation_initial_conditions`` already did. Only shock-free models
  are affected: when shocks exist the value is non-empty and the
  normalization does not fire, so every other path is byte-identical.

- **``solve_prune_hank_constants``: exact pruning of the variables that are
  constant in the linearized model** (new option, default ``false``; sideline
  development pending merge). A binary, parameter-free switch: the solver
  detects variables whose equations pin their first-order deviations at
  exactly zero at the evaluated steady state -- canonically the rail-stuck
  HANK interpolation weights (saturated clamps), but also slack borrowing
  multipliers, pinned policies, and fixed aggregates; the mechanism is
  generic -- removes them together with their defining equations by constant
  propagation on the evaluated Jacobians, solves the smaller system through
  the unchanged first-order machinery, and re-inflates the solution with
  exact-zero rows and columns. The answer is identical with the option on or
  off (measured IRF gaps 1e-13..1e-10 across the HANK replication suite);
  only the cost of solving it changes. The reduction is on the SYSTEM and
  happens ahead of whichever solver runs -- a direct decomposition for a
  constant-parameter model, the iterative MSRE solvers under regime
  switching -- so it is tied to no particular solution method. Under regime
  switching only variables trivial in every regime pair, at equal
  steady-state levels, are removed.
  If the reduction would not be square it refuses loudly and solves the full
  system. The option is HANK-named because that is where it is validated:
  the detector itself is model-agnostic, and it has been exercised across
  the model suite including models with no heterogeneity.
  (``+rise/+engine/+dsge_tools/+solve/{condense_ first_order,inflate_first_order}.m``, hook in ``first_order.m``.)
  Applies at first order only; higher orders continue to run at full
  dimension on the inflated ``Tz``.

- **Stability check: inert states dropped before the Kronecker build**
  (``+rise/+engine/+dsge_tools/+stability/is_stable_system.m``), gated on
  ``solve_prune_hank_constants``. States whose transition rows are exactly
  zero in every regime and solution are identically-zero directions;
  mean-square stability is invariant to removing them, and the Kronecker
  build is quadratic in the state count. This makes the MSS check feasible
  for pruned heterogeneous-agent models under regime switching: the
  2-regime Auclert HANK previously ran out of memory in
  ``gupta_murray_hassibi``; with ``solve_prune_hank_constants`` it now checks
  clean in seconds. Gated rather than unconditional so that every path not
  asking for pruning stays byte-identical -- without pruning an inert
  direction carries solver round-off (~1e-8) rather than an exact zero, so
  the scan would cost a pass and find nothing.

- **HANK: ``span`` declarations are now guarded, twice.** A ``span = s`` bracket
  window is an assertion about where the savings policy lands, and its
  validity depends on the parameterization -- so it is now checked
  machine-side instead of trusted. (1) The steady-state engine measures the
  landing distance ``d(i,j) = |b*(i,j) - i|`` at every solve on the full grid
  (``out.landing``: ``dmax``, ``min_safe_span``, per-cell ``d``), so the safe ``s``
  is handed to the user for free rather than guessed; if a declared span is
  violated at the current parameter values the engine stops with a hard
  error (``het:engine:spanViolated``) naming the violation and the minimal
  safe declaration -- the span-pruned model would be well-formed and could
  solve cleanly to a wrong steady state, so nothing downstream would catch
  it. (2) The parser now rejects ``span = 0`` (and non-integer or negative
  spans) with a clear message (``hank_expand:span``); previously it crashed
  with a raw index error, since brackets run ``1..N-1`` and the top cell's
  window would be empty. Validated in
  ``peer-feature-mining/HANK in Dynare/_experiments/exp29_span_guard.m``:
  span-free and ``span=1`` Auclert steady states identical to 1e-12; a
  uniform grid plus deep impatience (``beta=0.60``) trips the guard with
  "lands up to 3 brackets away (10 of 20 cells); minimal safe span=4".

- **HANK: generic steady-state engine**
  (``+rise/+engine/+dsge_tools/+sstate/+het/{engine,steady_state}.m``, new). A
  nested solver for the heterogeneous block, driven entirely by the parsed
  heterogeneity context: household policy by time iteration (all grid cells
  solved simultaneously by a vectorized Newton, Fischer-Burmeister at the
  borrowing constraint), cross-section by Young (2010) forward iteration,
  aggregation through the model's ``@Agg`` definitions, wrapped in a damped
  price/aggregate fixed point (``spec.prices`` closure) and an outer Broyden
  calibration (free parameters against named or function-handle targets;
  ``free = {}`` for pure fixed-point models a la Krusell-Smith). Validated
  against three independent hand-rolled steady states: Auclert et al.
  one-asset (all ten digits, full and span-pruned variants), McKay-
  Nakamura-Steinsson (1e-12 on every calibrated value), Krusell-Smith
  (imposed residual 5.6e-12 over the full expanded system), and
  Bayer-Luetticke (5,625 equations, borrowing to -100, entrepreneur via a
  bin indicator; imposed residual 1e-12, ~30x faster than the reference
  solver). Hardened for that last class: a guarded power (``het_pow``;
  negative exponents floor the base at 1e-12 -- the reference solvers'
  marginal-utility convention, inert on interior solutions) so "dead"
  grid cells with infeasibly negative consumption stay real; NaN-loud
  convergence tests (MATLAB's max ignores NaN, which had let a poisoned
  solve report success); per-cell monotone backtracking; and a Broyden
  that never commits a non-finite step near an asset-demand asymptote.

- **HANK parser: correctness fixes and new constructs** (``hank_expand.m``):

  1. *Individual leads on the left-hand side are now expanded correctly.*
     Previously a lead written on the LHS (e.g. ``beta*Va{t+1} = ...``) was
     silently renamed to the own-cell next-period variable -- no idiosyncratic
     expectation, no policy interpolation -- and solved a different model.
     Both sides now go through one expansion path; lead-free sides reduce to
     the old rename byte-for-byte.
  2. ``@Agg`` accepts nested parentheses (``@Agg(c{t}^(-sigma))``); the argument
     is scanned to its matching close paren. Leads/lags inside ``@Agg`` are
     rejected at parse time (previously mis-expanded silently).
  3. *Bin indicators* ``<axis>_at_<j>`` resolve at expansion time to literal
     0/1 by grid position -- addressing a bin by position where the bare
     coordinate only carries its value (per-type income, e.g. the
     Bayer-Luetticke entrepreneur). Out-of-range indices are parse errors.
  4. *Bracket window* ``@heterogeneity_axis(a, N, policy = ap, span = s)``
     limits the live interpolation weights to ``s`` brackets either side of
     each cell, cutting the expansion from O(N^2) to O(N) endogenous
     variables (measured: 7692 -> 1420 at N = 60). A window covering the
     grid reproduces the default expansion byte-identically.
  5. The parsed heterogeneity context (axes, cells, name maps, raw equation
     classification) is persisted to ``model_data.heterogeneity``; it is ``[]``
     for models without heterogeneity constructs.

- **DSMH merged and corrected: ``rise_dsmh`` + ``rise_dsmh2`` are now one
  sampler** (``classes/stats/mcmc/rise_dsmh.m``; ``rise_dsmh2`` is REMOVED --
  its option names carry over unchanged, so a caller just renames the
  function). The merged implementation keeps the usrsmplr-compatible
  interface and fixes/adds the following:

  1. *Striated-move acceptance is now the equi-energy ratio.* The old
     ``rise_dsmh2`` accepted a resample from the previous stage's energy
     band with the plain tempered ratio ``exp(lambda_i * dlogp)``, although
     the proposal is drawn from the PREVIOUS tempered target restricted
     to the shared striation; the correct Hastings ratio is
     ``exp((lambda_i - lambda_{i-1}) * dell)`` (Kou-Zhou-Wong 2006 /
     Waggoner-Wu-Zha 2014). Wrong occupancy shares on multimodal
     targets under the old ratio. Striated jumps now also recycle the
     stored posterior values, so they cost zero new evaluations.
  2. *Out-of-bounds proposals are rejected, never clamped* (same defect
     fixed for rwmh/imh/apt on 2026-07-05): with ``lb``/``ub`` supplied, the
     raw proposal is tested BEFORE the frame is applied, which is exact
     MH on the bound-truncated target even when a legacy clamping frame
     is passed in.
  3. *Covariance adaptation stores the lower Cholesky factor.* Both old
     files assigned the upper ``chol`` factor to the variable used as a
     lower factor, so the adapted proposal had the right spectrum but a
     rotated orientation (efficiency loss only -- the proposal stays
     symmetric, so no bias; ``rise_dram`` had the same quirk and is fixed
     too). Between stages the scalar scale now also adapts toward a
     0.234 acceptance target.
  4. *Likelihood tempering* (``logprior`` option): the ladder tempers the
     likelihood only and the prior enters every stage untempered -- the
     Waggoner-Wu-Zha construction; the default remains posterior
     tempering. The default ladder is now geometric ``(i/H)^psi``
     (psi = 4, dense near lambda = 0) instead of uniform.
  5. *Stepping-stone marginal-data-density by-product*
     (``extras.log_mdd``): with a lambda = 0 anchor (iid uniform-box
     draws given ``lb``/``ub`` under posterior tempering, or iid prior
     draws via ``prior_sampler`` under likelihood tempering) the ladder's
     normalizing-constant increments are accumulated into a log-MDD
     estimate at no extra sampling cost beyond the anchor stage. Each
     increment is estimated two-sided (forward from the lower rung's
     draws, reverse from the upper rung's, averaged in the log domain),
     which cancels the leading-order Jensen bias of one-sided
     stepping-stone estimation.
  6. *Exact checkpoint/resume* (``stop_after_sweeps`` / ``state`` options):
     the full sampler state including the RNG is returned and can be
     passed back; a resumed run reproduces the uninterrupted run draw
     for draw. ``burnin``, previously ignored by both implementations, is
     now a per-stage burn-in.

  All of it is locked by the new suite
  ``RISE-unit-tests/infrastructure/samplers/mcmc/DsmhCorrectnessTest.m``
  (analytic-moment targets, a mode-on-the-bound target with a clamping
  frame, a 0.7/0.3 bimodal mixture started in the minor mode with a
  plain-MH control, analytic-evidence checks for both tempering modes,
  and a draw-for-draw resume test).

- **Fix: ``rise_dram`` adapted proposal used the upper Cholesky factor**
  (same orientation quirk as item 3 above; efficiency, not bias). The
  adaptation now stores the lower factor, and the delayed-rejection
  correction ``iL`` is consistent with it.

- **Performance: ``rise_dates.dates`` arithmetic vectorized** (``@dates/dateshift``,
  ``@dates/double``). ``dateshift`` (called by ``plus``/``minus``/``colon``) and ``double``
  looped element-by-element, invoking MATLAB's already-vectorized ``datetime``
  ``dateshift``/``datenum`` once *per element*. They now batch the underlying
  ``datetime`` call over the whole array. Results are bit-identical across all
  frequencies (D/W/M/Q/H/Y; char, serial, ``+``, ``-`` all verified against the prior
  implementation) and the 234-test ``infrastructure/time_series`` area stays fully
  green. Measured ~3.6x faster on a monthly range build + ``double`` and ~2.7x on a
  quarterly ``colon``. First cut of a broader effort to keep dates/time-series off
  the estimation hot path (the per-evaluation likelihood loop is already
  double/struct-only; the remaining date cost is one-time data assembly and
  filter-output serialization).

- **Fix: stochastic regime draws in simulation restored** (regression introduced
  2026-06-12 by ``45255f15``, which made the OBC constrained likelihood deterministic
  by selecting the modal next regime in ``compute_transition`` -- correct for
  reproducible mode-finding, but the routine is shared with plain simulation via
  ``forecast_unconstrained``, so ``simulate()`` paths silently froze in their initial
  regime whenever chains are persistent). Regime transitions in simulation are
  drawn from the transition matrix again. Design: stochastic is now the DEFAULT
  semantics of ``compute_transition``; the constrained/violation-check forecasts
  (``forecast_with_constraints``, ``forecast_with_constraints_functional``) explicitly
  opt out (``stochastic=false``) to keep the OBC likelihood bit-reproducible
  (verified). New option ``simul_stochastic_regimes`` (default ``true``) lets a user
  reproduce the deterministic modal path in simulation if desired. Note: GIRFs
  again average over stochastic regime futures, the textbook definition.

- **Fix: numeric ``simul_regime`` honored again.** The documented scalar/vector
  forms of ``simul_regime`` (impose a regime path in simulation) were silently
  dropped (``simul.load_options`` only handled the function-handle/switch-rule
  forms; ``set_simulation_regimes`` had become orphaned). A numeric path is now
  expanded in ``set_simulation_initial_conditions`` and overrides the regime draw;
  imposed paths are reproduced exactly.

- **Exact (Durbin-Koopman) diffuse initialization** for the Kalman filters, selected by
  ``kf_init_variance='exact'`` (the same option that takes ``[]`` for the unconditional/Lyapunov
  default or a positive scalar for Harvey big-kappa -- exact is a third value there, not a
  separate option). The diffuse directions are resolved in
  a short exact diffuse phase, so the likelihood is kappa-free -- it drops the spurious
  ``-q/2*log(2*pi*kappa)`` offset that big-kappa carries and that contaminates MLE, marginal
  data density and model comparison. ``initialization.m`` builds ``P_inf`` on the nonstationary
  (unit-root) subspace of the state transition; ``imm.m`` runs a guarded two-covariance
  diffuse phase; single-regime ``exact`` models are routed to ``imm``, so it works for
  constant-parameter and regime-switching models alike, and pinning a non-``imm`` filter with
  ``exact`` warns. The default big-kappa path is unchanged (filtering 47/47, obc 23/23).
  Validated ``ll_exact = ll_bigk + q/2 log(kappa)`` at O(1/kappa). See the "Diffuse
  initialization" docs page and the ``Estimation/diffuse_initialization`` tutorial. Exact
  diffuse is available in ``imm`` (default) and in the GPB-family filters ``immn``,
  ``kim_nelson`` and ``gpbn`` (pin via ``kf_user_algo``); all four give the same diffuse loglik.
  Pinning a filter that does not implement it (``imm_o``, the constrained cell filter) warns
  and falls back to big-kappa. The diffuse subspace is now the nonstationary **invariant
  subspace** (ordered real Schur), correct for **multiple and defective** unit roots -- e.g.
  a local linear trend (``[1 1;0 1]`` Jordan block, 2-D diffuse subspace, one eigenvector) --
  not just a random walk; ``P_inf`` is carried through the filter's state-space augmentation.
  Validated across all four filters on random-walk (q=1) and local-linear-trend (q=2)
  models (``testDiffuseInit``, 12 cases); default big-kappa path unchanged (GPB suite 13/13).
  Two guards keep exact diffuse within its valid scope (the regular, regime-invariant
  diffuse block): if the nonstationary invariant subspace differs across regimes
  (regime-DEPENDENT diffuse block), or if measurement noise ``H`` is non-diagonal (the
  univariate diffuse cascade assumes diagonal ``H``), ``initialization.m`` warns and falls
  back to big-kappa rather than compute a silently-wrong likelihood.
  Exact-diffuse **smoothing** is wired into **all four** filters (imm, immn, kim_nelson,
  gpbn): the forward pass stores the predicted ``P_inf`` for the diffuse steps and the
  smoother (``utils.filtering.gpbn.smoother`` for imm/immn/gpbn, ``utils.filtering.smoother``
  for kim_nelson) runs a Durbin-Koopman **dual-accumulator** (``r0``,``r1``) backward
  recursion over the diffuse span, so smoothed states over the first ``d`` periods are
  exact (previously they were degraded). Guarded: the branch fires only when ``P_inf``
  storage is present, so every non-diffuse smoother is byte-identical. ``simple_kalman_step``
  now returns the predicted ``P_inf`` and the real diffuse-phase innovation. Validated:
  smoothed states match big-kappa at O(1/kappa) for all four filters (local linear trend,
  q=2); ``testDiffuseInit`` 16/16.

- **The default regime-switching filter is now the Interacting Multiple Model
  (IMM) algorithm** (Blom-Bar-Shalom), replacing
  ``constrained_regime_switching_kalman_filter_cell`` for order-1 switching models.
  IMM's mixing is more general than the Kim-Nelson collapse the old default used.
  To reach parity, the two features the old filter had and ``imm`` lacked were
  ported into ``imm`` via new shared helpers so any switching filter can reuse
  them: ``rise.engine.filtering.filter_exo_driver`` (the observed-exogenous
  time-varying-transition driver, threaded into
  ``evaluate_filter_transition_matrix``) and ``rise.engine.filtering.filter_cond_data``
  (the real-time / observed-expectations conditioning slice, passed as the
  predict closure's 5th argument). ``imm`` also now reads ``data_y`` with future
  pages and the real-time ``syst`` fields. Filtering, OBC, regime-switching,
  estimation (MLE) and dynare-comparison suites are all green under the new
  default. The GPB-family filters (``gpbn``, ``imm_o``, ``immn``, ``kim_nelson``) can
  adopt the same two helpers to gain these features. The previous filter remains
  available via ``kf_user_algo``.

- **Leads greater than one are now reduced without violating Jensen's
  inequality -- everywhere, by default.** The parser's lead reduction previously
  chose between two mechanisms via the ``violate_jensen_inequality`` option: the
  default atomic chain (fast, but ``E[g(x)] -> g(E[x])`` for nonlinear leads) or a
  whole-equation backshift that closed each reduced equation with a martingale
  row ``AUX{t+1}=0``. The backshift was Jensen-correct but its martingale row has
  no contemporaneous variable, making ``A0`` singular and breaking iterative
  solvers (cyclic reduction; the ``mfi`` optimal-policy iteration) -- so it could
  never be the default. Both are replaced by a single term-scoped reduction
  (``get_rid_of_excess_leads_preserving_jensen``): a term that is NONLINEAR in a
  forward variable is wrapped in an auxiliary ``z{t}=shift(term,-1)`` and replaced
  by ``z{+1}`` (carrying the whole term's one-period-ahead conditional
  expectation, ``E[g(x)]``), while AFFINE leads use the ordinary atomic chain
  (``chain_excess_leads``). Every auxiliary keeps a contemporaneous variable, so
  ``A0`` stays non-singular. A recursive-descent classifier decides linearity;
  a lead scaled by a quantity known at ``t`` stays affine, a lead inside a
  function/power or multiplied by another forward is wrapped. Models with no
  nonlinear leads (all linear models, and derived optimal-policy FOCs such as
  fm95) are reduced byte-for-byte as before. The ``violate_jensen_inequality``
  option and the backshift reduction are **removed**. Full suite green (1004
  tests), including the ten lead/lag regressions the backshift used to break.
  Details and the closed-form ``sigma^2`` derivation:
  ``development/jensen_memo/jensen_inequality.pdf``.

- **New ``rise.microfound`` -- write the economics, not the algebra.** A model
  *definition* builder: declare the atoms and each agent's optimization problem
  (period objective, choice variables, named constraints) in MATLAB, and it
  derives the first-order conditions symbolically and produces a plain RISE
  model -- no ``@optimization_problem`` block, nothing hidden. The derivation
  reuses the FOC engine's beta-weighted differentiate-then-shift assembly via
  new shared primitives ``rise.engine.foc.{change_symbols,symbolic_differentiate, attach_lagrange_multipliers}`` (byte-for-byte extracted from the optimal-policy
  engine). Hand the builder to ``dsge_model`` three ways: ``dsge_model(p)``,
  ``dsge_model(code(p))``, or write a ``.rs`` with ``rs(p,file)``. microfound builds
  the DEFINITION only -- calibration (``solve(m,'parameters',struct)``) and
  solving are ``dsge_model``'s responsibility. Derived multipliers are named
  ``<agent>_lambda_<constraint>``. Switching parameters in an agent's problem
  yield correct probability-weighted switching FOCs automatically. IMPORTANT:
  the steady state is NOT provided analytically (that would mean hand-deriving
  what the tool automates) -- ``dsge_model`` forms the steady-state system from
  the derived FOCs and solves it numerically; ``steady_state(...)`` supplies only
  a rough initial guess. Validated by re-authoring 9 gEcon example models (rbc,
  rbc_hf, rbc_cu, rbc_ts, rbc_ic, home_production, rbc_mc with household-SDF
  discounting, tc_rs two-country risk sharing, ttb time-to-build) -- all derive
  and solve, and rbc_ts matches gEcon's own derived steady state. See
  ``examples/microfound/``. A **companion pre-pass** auto-reduces lags > 1 before
  the FOCs (a control at lag k becomes a chain of co-state auxiliary controls
  with identity constraints; a non-control variable becomes auxiliary endogenous
  with model identities) so no manual auxiliary chain is needed -- e.g. a
  2-period time-to-build ``K = (1-delta)*K{-1} + I{-2}`` derives and solves
  directly. Leads >= 2 in an agent constraint are reduced by a Jensen-preserving
  subexpression aux: an additive term ``T`` at lead >= 2 becomes ``z{+1}`` with the
  identity ``z = shift(T,-1)``, so ``z`` carries the WHOLE (possibly nonlinear)
  term's one-period-ahead conditional expectation (``E[g(x)]``, not ``g(E[x])``) and
  the constraint keeps its date (its multiplier coupling intact). Verified to
  reproduce the direct FOC ``(1/beta^2)*lambda_{t-2}*(dg/dx)`` exactly, and a
  nonlinear lead-2 model at second order matches the analytic Jensen gap
  ``exp(0.5*sig^2)`` (unit test ``models/dsge/microfound/testMicrofoundJensen``). A
  time-varying (stochastic) discount factor is handled correctly: the backward
  multiplier weight uses the LAGGED one-period discount ``1/SDF_{t-1,t}`` (for a
  constant discount this equals ``1/beta``, so the planner path is unchanged). HANK
  heterogeneity, Kuhn-Tucker inequality constraints, recursive (Epstein-Zin)
  preferences, and gEcon-style index-set templates are planned extensions.

- **``fix_point_iterator`` now speaks plain MATLAB option names internally --
  the last ``fix_point_*`` names are gone.** The shared fixed-point engine used
  to rename incoming options to internal ``fix_point_TolFun`` /
  ``fix_point_maxiter`` / ``fix_point_verbose`` / ``fix_point_explosion_limit`` /
  ``fix_point_valid_func`` fields. It now reads plain names directly:
  ``TolFun``, ``MaxIter``, ``Display`` (``'iter'``|``'none'``|``'final'``, with the
  optimset ``'off'`` and a logical still accepted for back-compat),
  ``ObjectiveLimit``, ``ValidFunc``. Behaviour at every call site is unchanged
  (the options structs already carried the plain names); the internal
  ``update_fix_point_options`` rename step is removed. Callers that constructed
  ``fix_point_*`` fields (the Maih-Waggoner ``update_options``, the Lyapunov
  ``doublingFixIter``) and the solver unit tests were updated to the plain
  names. Verified: solvers_generic 195/195, steady-state 262/262, and the
  Kalman/Lyapunov path (filter) unchanged.

- **New MSRE solver ``mafia`` -- Multisecant-Accelerated Functional Iteration
  Algorithm.** The Picard iteration behind ``mfi`` converges only on a
  contractive fixed-point map; on a genuinely non-contractive map (one
  Jacobian eigenvalue > 1) it diverges, returning retcode 21, even where a
  stable solution exists and ``mnk``/``mn``/``fwz`` find it. ``mafia``
  (``solve(m,'solver','mafia')``) iterates the SAME map as ``mfi`` -- so it keeps
  ``mfi``'s tendency to land on the stable minimum-state-variable solution
  rather than an explosive one -- but accelerates it with self-tuning windowed
  Anderson acceleration (a multisecant/quasi-Newton scheme using only
  evaluations of the map, no Jacobian, with QR condition-based regularisation
  and a non-monotone-bounded safeguard). It has NO tuning knobs: one fixed
  configuration converges the easy-contracting, non-contractive and
  slow/ill-conditioned regimes alike, so it can be pointed at estimation and
  adapt per parameter draw without the caller knowing which regime each draw
  is in. On ``fin_frics_switch`` (Picard-map spectral radius 1.0163) plain ``mfi``
  returns retcode 21 while ``mafia`` recovers the ``mnk`` solution to ~1e-7, and it
  is faster than ``mfi`` on the maps ``mfi`` already solves. ``mfi`` and
  ``fix_point_iterator`` are unchanged (byte-for-byte). Theory and the
  parameter-free three-regime validation:
  ``development/MFI_Anderson_Acceleration_Note.pdf``,
  ``development/ADAPTIVE_AFI_VALIDATION.md``; regression:
  ``RISE-unit-tests/.../solvers_generic/testMafiaSolver.m``.

- **New option ``sstate_warm_start`` (default ``false``): warm-start the
  numerical steady-state solve from the last successful solution.**
  Profiling shows the numerical steady-state search
  (``rise.engine.rstst.solve_blocks``) is 42-69% of every model re-solve,
  and estimation re-solves the steady state from scratch on every draw
  even though successive draws are ~0.1% apart in parameter space. With
  ``sstate_warm_start=true``, the search is seeded with the last
  successful steady state of the same model (persistent cache keyed by
  a stable model fingerprint — the object itself is re-created per
  draw); the current draw's parameters always prevail, and a failed
  warm-started solve is automatically retried once from the cold
  default, so a stale cache entry can never make a solvable draw fail.
  With a unique steady state the warm and cold solves land on the same
  solution (verified to <1e-8 in the tests). Measured on 20 jittered
  re-solves (MCMC-like): fs2000 solves 2.2x faster end-to-end
  (order 2); models whose cold default already sits at the steady state
  (e.g. linear switching models with a zero steady state) see no
  change. Reset the cache with
  ``rise.engine.rstst.warm_start_cache('clear')``. Default OFF: behavior
  is unchanged unless the option is set.

- **Retired the global ``fix_point_*`` options** (``fix_point_TolFun``,
  ``fix_point_maxiter``, ``fix_point_verbose``, ``fix_point_explosion_limit``,
  ``fix_point_valid_func``). These were a
  single shared knob that silently affected every consumer of the internal
  fixed-point iterator (MSRE solvers, optimal policy, Lyapunov doubling,
  linear-systems solves, the mdd bridge sampler, ...). Each call site now
  reads its own scoped option with unchanged defaults: the solver options
  cell (``solver={'mfi','MaxIter',...,'TolFun',...}``; ``Display='iter'``
  replaces ``fix_point_verbose``), ``sstate_solver{2}``,
  ``solve_linsyst_solver{2}``, ``lyapunov_algo{2}``, and
  ``mdd_bridge_TolFun``/``mdd_bridge_maxiter`` for the bridge sampler. The
  discretion iteration cap is ``solve_discretion_maxiter``. Passing a
  ``fix_point_*`` option now errors as an unknown option.

- **MCMC correctness fixes (three defects) + a procedural-correctness
  test suite.** An audit of the sampling machinery for algorithmic
  correctness (acceptance ratios, proposal corrections, bounds handling,
  sign conventions) found and fixed three defects:

  1. *Independence MH drew from the wrong distribution.* ``one_step_mh``
     (and the delayed-acceptance variant) proposed from a random walk
     centered at the current state while the acceptance ratio applied the
     independence-sampler density correction with the fixed proposal
     mean — the chain's stationary distribution was ``pi(x)/q(x)`` instead
     of ``pi(x)`` (inflated tails). IMH now draws from the fixed proposal
     distribution, as the correction assumes.
  2. *Out-of-bounds proposals were clamped, not rejected.* RWMH/IMH/APT
     proposals were clamped onto ``[lb, ub]`` by the recentering frame and
     evaluated there, while the acceptance ratio assumed a symmetric,
     uncorrected proposal — piling probability mass on active bounds.
     Proposals now use a constraints-only frame (``reframe`` style 0) and
     out-of-bounds draws are rejected in the accept/reject step, which is
     exact MH on the bound-truncated target. Behavior is unchanged when
     the bounds do not bind (identical RNG stream, no redraws).
  3. *DRAM's delayed-rejection stage used a mis-scaled proposal density.*
     In ``rise_dram``, the second-stage acceptance's ``q_ratio`` omitted the
     ``1/adascale^2`` scaling of the first-stage proposal covariance, so
     the DR acceptance probability was wrong whenever ``adascale ~= 1``
     (the default is ``2.38/sqrt(d)``).

  All three are locked by the new suite
  ``RISE-unit-tests/infrastructure/samplers/mcmc/McmcCorrectnessTest.m``,
  which checks each sampler (RWMH, IMH, APT, slice, DRAM with the DR
  branch actually enabled, shortcut wrappers) against targets with
  analytic moments, including a target whose mode sits ON a bound.

- **New option: ``sstate_file_is_authoritative`` (default ``false``).** When a
  user-supplied steady-state file (or steady-state model) returns a nonzero
  retcode — "no solution for these parameters" — the engine used to ignore
  the verdict and fall through to the full numerical steady-state search,
  a 1-10 s tax per bad draw during MCMC estimation. With the option set to
  ``true``, the file's verdict is trusted: its retcode propagates immediately
  and the numerical search is skipped. A large residual with a zero file
  retcode still triggers the numerical search as before. On the
  ``sstate_loop`` path the verdict is checked at the initial loop point only.
  Default ``false`` preserves the previous behavior exactly.

- **New: ``texreport(m)`` — paper-ready LaTeX model appendix (headliner 4,
  slice 2).** One command writes a compilable ``.tex`` (or an ``\input``-ready
  fragment): the model equations rendered as mathematics by a new
  RISE-to-LaTeX translator (``rise.engine.reporting.eq2tex`` — a
  recursive-descent parser emitting ``\frac`` for divisions, braced powers,
  operator names for known functions, time subscripts ``x_{t+1}`` from the
  ``{t+1}``/``{+1}`` forms, ``\bar{x}`` for ``{stst}``, and greek symbols via the
  model's tex names or a built-in table), plus calibration and
  steady-state tables per regime (booktabs). Auxiliary LEAD\_/LAG\_
  equations omitted by default; translation failures fall back to
  verbatim rather than dying; multiple parameterizations are sliced into
  one appendix each; ``'compile', true`` runs pdflatex. Unit-tested down to
  exact translator strings, all 16 fs2000 equations translate with
  balanced groups, and the compile test produces a real PDF.

- **New: ``webreport(m)`` — a self-contained interactive HTML model report.**
  One command produces ONE file — no external assets, no server, opens
  offline in any modern browser, shareable with coauthors who do not have
  MATLAB. Tabs: overview (regimes, solver, root moduli for
  constant-parameter models), the model equations, calibration and steady
  state per regime, an **interactive IRF explorer** (shock selector,
  variable picker, one panel per variable with the regimes as series,
  crosshair hover reading out every regime at the hovered horizon, regime
  legend and end-of-line labels), and a data-table tab carrying every
  plotted number. Light and dark mode follow the system preference;
  charts use a colorblind-validated palette. Engine:
  ``@rise_model/webreport.m`` + the template in
  ``+rise/+engine/+reporting/webreport_template.html``; numbers are
  embedded as JSON inside the file. First slice of wahoo headliner 4 —
  posterior/estimation panels and a LaTeX model-appendix export are the
  planned follow-ups.

- **Failure diagnosis phase D5 — the return-code registry now speaks.**
  The terse registry messages behind ``decipher`` were rewritten to state
  the cause and the next step: ``'System unstable'`` became "No
  mean-square-stable solution: the switching system's second moments
  explode under the regime process (this is the MSS test, not an
  eigenvalue count) — run diagnose(m) for the MSS radius and the culprit
  regime"; ``'Explosive solution'``, ``'Nans in solution or no solution'``,
  the steady-state and simulation codes similarly. Codes, names and
  categories are unchanged — only the human-facing text. A *Failure
  diagnosis* chapter joins the docs (Working with a model) and a runnable
  tutorial joins rise-modern-tutorials (three deliberately broken models,
  the homotopy stall-as-diagnosis, and the MSS-vs-conditional-explosiveness
  lesson).

- **New: failure diagnosis phase D4b — perfect-foresight shock-size
  continuation.** A new ``simul_homotopy`` option (default false): when the
  stacked perfect-foresight solve fails at the full shock size, the
  exogenous path is walked up from zero (adaptive step, each trial
  warm-starting the free endogenous cells from the previous step's path;
  initial/terminal conditions untouched), so the full-size problem is
  approached from a nearby solution instead of the steady-state guess it
  just failed from. If even the ladder stalls, the reached share of the
  shock size is itself the diagnosis ("solvable up to 60% of the shock —
  the full-size shock likely pushes the path where no bounded solution
  exists"), surfaced through ``diagnostics_on_failure``. Inner ladder trials
  keep their per-step diagnostics silent. Also fixes an edge case in the
  D3 worst-cell analyzer when every period's residual is non-finite.
  Honest test note: the mechanism is verified by traced runs and a
  no-interference test; a genuine cold-failure integration fixture is
  still missing because the stacked Newton solves the available smooth
  models even at absurd shock sizes.

- **New: failure diagnosis phase D2 — the stability analyzer, with the
  right vocabulary per model class.** When the solve fails, ``diagnose(m)``
  now explains why in the language that applies: for CONSTANT-PARAMETER
  models, Blanchard-Kahn root counting — "6 stable roots but only 5
  predetermined variables: indeterminacy", with the root table (the failed
  solve now retains its eigenvalues for exactly this purpose) and a
  clearly-labeled Taylor-principle heuristic; for REGIME-SWITCHING models,
  MEAN-SQUARE STABILITY only — the MSS spectral radius vs the criterion,
  per-regime conditional radii and persistence, and the culprit regime
  ("regime 2 is conditionally explosive and too persistent") — never an
  eigenvalue count. Iterative MSRE solver failures (no convergence / NaNs /
  divergence) get their own narrative with solver alternatives. Engine:
  ``rise.engine.diagnose.stability``; ``model_data`` gains a template-declared
  ``diagnostics`` field (keeps struct arrays homogeneous under estimation).

- **New: failure diagnosis phases D3+D4a — no more silent NaNs, and a
  homotopy fix-action.**

  - *The silent-NaN era ends.* A diverging simulation used to return a
    NaN-filled database without a word. Now the failure site reports the
    first non-finite value — period, variable(s), regime — the growth
    profile of the run-up, and the pruning hint when running above order 1
    unpruned (``rise.engine.diagnose.simulation``). A new
    ``diagnostics_on_failure`` option ({'warn'} | 'error' | 'off') governs
    surfacing everywhere; 'off' restores the legacy silence.
  - *Perfect-foresight worst cells.* When the stacked solver fails to
    converge, the stacked residual at its final iterate is reshaped to
    (equation x period) and the worst violations are reported, plus the
    periods carrying 80% of the residual mass — localizing the failure to,
    say, "equation 3 around periods 12–15" instead of "did not converge"
    (``rise.engine.diagnose.perfect_foresight_stack``; the existing
    ``multiplicity_diagnostic`` is suggested for uniqueness questions).
  - *Homotopy fix-action.* Every successful steady-state solve now
    remembers its calibration; ``[report, mfixed] = diagnose(m, 'fix', "homotopy")`` walks the parameters gradually from that last working
    calibration to the failing target, warm-starting each step from the
    previous solution (with a cold fallback when a warm start is unusable,
    e.g. a log variable whose steady state went non-positive). On success
    the model comes back solved at the target; on stall the breakpoint is
    itself the diagnosis — "solvable up to 47% of the way: the target
    calibration likely admits no steady state" — and the model is returned
    solved at the breakpoint mix for inspection.

- **New: ``diagnose(m)`` — failure diagnosis, phase 1 (steady state).** The
  first slice of the failure-diagnosis suite
  (``development/FAILURE_DIAGNOSIS_MEMO.md``): instead of a bare
  ``Steady state could not solve``, ``diagnose(m)`` reports the worst
  equations (per-equation residuals against the engine's 1e-6 acceptance
  floor, labeled by equation tag or the equation text itself), the
  variables and parameters entering them, a structural-rank check of the
  steady-state incidence (a declared-but-never-used variable is named
  outright), and concrete next steps. Diffuse failures (most equations
  above the floor) are called out as such — one global culprit is more
  likely than the top-ranked equations — rather than overattributed.
  When the steady state is fine but the solve fails, the deciphered
  return code is surfaced with model-class-aware stability language
  (eigenvalue/Blanchard-Kahn reasoning only for constant-parameter
  models; mean-square stability for regime-switching ones); the full
  stability analyzer is the next phase. Engine:
  ``rise.engine.diagnose.steady_state``.

- **Fix: ``load_solution`` errored when called with a model object.**
  ``rise.engine.dsge_tools.solve.load_solution`` documents that it accepts a
  model object or a ``model_data`` struct, but it read ``obj.topology`` (and,
  further down, ``obj.exogenous``) — fields the class facade does not
  forward as dependent properties — so a direct call on a ``dsge_model``
  errored. The function now normalizes a model object to its ``model_data``
  after the solve step; both calling conventions are regression-tested and
  deliver identical solutions.

- **New: unified ``sensitivity()`` front-end — one design, three GSA views.**
  A ``sensitivity`` method on all model classes (engine:
  ``rise.engine.surrogates.sensitivity_report``) runs the whole
  global-sensitivity toolbox on ONE evaluated design: analytic PCE-Sobol'
  main/total indices (who drives the variance of the log posterior over
  the prior), HDMR first-order indices from an independently fitted
  metamodel on the same design (a built-in cross-check, validated against
  analytic benchmarks), and Monte-Carlo filtering (who decides whether a
  draw lands in the top region of the target, measured by per-parameter
  Kolmogorov-Smirnov separation with asymptotic p-values — the full ``mcf``
  object is returned for its plots). Everything lands in one table;
  HDMR/MCF recycle the emulation design rather than sampling again. Note:
  only HDMR's first-order column is reported — the legacy engine's
  level-2+ aggregates were found to carry systematic leakage/overcounting
  at realistic sample sizes (e.g. reading 0.24 on a null pair, 1.19 on a
  pure pair with share 1.0); investigating that estimator is left as a
  follow-up, and interaction/total readings come from the PCE columns.

- **New: adaptive surrogate tools — refinement, EGO mode finding, active
  subspaces.** Three tools close the surrogate lifecycle on the same
  ``fit``/``predict`` contract. ``refine(s, f, nadd)`` adds true evaluations
  where the emulator is least sure (the GP's pointwise sd; space-filling
  for the PCE) and refits in rounds. ``rise.engine.surrogates.ei_maximize``
  runs the classic EGO loop (Jones-Schonlau-Welch expected improvement) to
  locate the optimum of an expensive function — e.g. a posterior mode — in
  a few dozen true evaluations (finds the Branin optimum within a
  60-evaluation budget in the test suite), returning the fitted GP for
  reuse as a delayed-acceptance screen. ``active_subspace(s)`` estimates the
  gradient-covariance eigendecomposition under the prior with gradients
  taken on the surrogate (zero additional true evaluations), exposing the
  few linear parameter combinations the posterior actually responds to,
  plus per-parameter activity scores.

- **New: Gaussian-process (Kriging) emulator with calibrated pointwise
  uncertainty.** ``emulate(m, 'method', "gp")`` — the second emulator on the
  same ``fit``/``predict`` surrogate contract. Detrended (universal) Kriging: a
  polynomial mean (``trend`` option; ``emulate`` uses ``"quadratic"``, the natural
  shape of a log posterior) plus an anisotropic squared-exponential GP on
  the residuals, inputs standardized by the design data (robust to
  fat-tailed priors whose mass is a sliver of the bound box),
  hyperparameters by marginal likelihood with restarts, analytic
  leave-one-out ``Q2`` diagnostics. Unlike the PCE's single global
  uncertainty number, ``predict``'s second output is **pointwise**: small
  near design points, growing away from them — the ingredient for adaptive
  designs and surrogate-guided optimization. On ``fs2000``: ``Q2`` ~0.92 from
  135 design points (order-2 PCE: ~0.81 from 165). No analytic Sobol'
  (use the PCE for sensitivity); both emulators work as delayed-acceptance
  screens.

- **New: chain checkpointing and resume for the MH samplers.** ``rwmh`` and
  ``imh`` accept ``store_file`` / ``store_every`` / ``resume``: the complete chain
  state (draws, current position, rng state, adaptive-tuning state, and
  the delayed-acceptance surrogate value when active) is checkpointed
  atomically to MAT v7.3 (HDF5), so long runs survive crashes. Resuming
  restores the rng and tuning states exactly, making the continued chain
  **bit-for-bit** identical to the uninterrupted one (unit-tested); a
  completed run can be extended by resuming with a larger ``N``. With
  ``nchain > 1`` a per-chain suffix is appended. ``apt``/``slice`` reject the
  option explicitly.

- **New: delayed-acceptance MCMC — exact posterior at a fraction of the
  true-likelihood cost.** The Metropolis-Hastings samplers (``rwmh``, ``imh``)
  gained a ``da_surrogate`` property: every proposal is first screened with a
  cheap surrogate of the log posterior (e.g. ``@(x)predict(s,x)`` from
  ``emulate``); only survivors pay for a true posterior evaluation, and a
  second accept/reject step (Christen & Fox, 2005) keeps the chain's
  stationary distribution **exactly** the true posterior — a poor surrogate
  costs acceptance rate, never correctness (unit-tested with a deliberately
  biased surrogate). Composes with the existing adaptive tuning and
  constraint machinery, and ``stats.funevals`` reports the true evaluations
  actually paid for. On ``fs2000``, ~70% of iterations are settled by the
  screen and wall time drops ~3x with indistinguishable posteriors. A
  surrogate that fails at the chain's initial point disables the screen for
  that chain (with a warning) instead of freezing it; the ``apt`` and ``slice``
  samplers reject the option explicitly. Example:
  ``rise-modern-tutorials/Estimation/surrogate_accelerated_mcmc/howto.m``;
  docs: the *Surrogates* page. Phase 1b of the
  surrogate-accelerated-estimation roadmap.

- **New: surrogate (emulator) layer — polynomial chaos with analytic Sobol'
  indices.** A new ``rise.engine.surrogates`` package (``pce``, ``doe``, abstract
  ``surrogate`` contract) plus an ``emulate`` method on all model classes. One
  call — ``[s, design] = emulate(m)`` — fits a polynomial-chaos expansion of the
  log posterior kernel (or of any user function of the parameters) from a few
  hundred true evaluations; afterwards ``predict(s, theta)`` approximates it in
  microseconds (~1500x per evaluation on ``fs2000``). The polynomial basis is
  derived **per parameter from its prior** (Hermite for normal priors,
  cdf-transformed Legendre otherwise), the design of experiments is drawn from
  the priors' inverse cdfs, and — because the basis is orthonormal under the
  prior — ``sobol(s)`` returns the exact variance decomposition of the emulated
  quantity (main and total Sobol' indices) analytically from the fitted
  coefficients, with no additional sampling. Leave-one-out cross-validated
  ``Q2`` and hold-out validation helpers included. Example:
  ``rise-modern-tutorials/GSAandUQ/UQ_surrogates/howto.m``; docs: the
  *Surrogates* page of the GSA/UQ chapter. This is Phase 1a of the surrogate-accelerated-estimation roadmap
  (``development/WAHOO_ROADMAP.md``); delayed-acceptance MCMC and a
  Gaussian-process emulator on the same contract come next.

- **Perfect-foresight solver: faster and more robust.** The deterministic
  (``perfect_foresight``) solver gained several improvements, all
  behaviour-preserving at their defaults:

  - the default ``sparse`` solver (``rise_newton``) now judges convergence on the
    **max-abs (infinity-norm)** residual instead of the 2-norm. The 2-norm of a
    stacked residual grows like ``sqrt(periods*endo_nbr)``, so a fixed tolerance
    was unreachable for long horizons even when every equation residual was
    tiny; the max-abs criterion is per-equation and size-independent.
  - a new **``JacobianUpdate``** option on the ``sparse`` solver enables
    modified-Newton **factorization reuse**: factorize the stacked Jacobian
    once and reuse the LU factors (``inf``), or refresh every ``K`` iterations
    (``K>1``); ``1`` (default) is full Newton. Reuse refreshes adaptively when the
    line search backtracks, and retries with a fresh factorization on a stall.
    On transition problems this leaves the solution unchanged while
    substantially cutting the time spent in the linear solve.
  - the per-period Jacobian block is now assembled directly as sparse
    (dropping a dense scatter + ``find``), and the line search no longer rebuilds
    the Jacobian on every backtracking trial.
  - a new **``simul_reuse_sstate``** option (default ``false``) skips the internal
    steady-state re-solve when the model is already solved, for pre-solved /
    repeated / extended-path runs.

- **Steady-state solve: false "could not solve" failures fixed.** The
  steady-state success verdict compared the residual against the solver's
  iterate tolerance (``sqrt(eps)``), which sits right at the achievable residual
  floor — so well-converged solutions on some models (e.g. the canonical
  ``ramst``) were falsely reported as failures. The acceptance verdict is now
  decoupled from the solver tolerance via a floor (``1e-6``, still tighter than
  typical practice and well above the floor), judged on the max-abs residual;
  a deliberately looser solver tolerance is still honoured, and achieved
  precision is unchanged. Separately, the engine now refreshes the parameter /
  definition slots from the current values on every solve, so parameter
  changes (estimation, calibration) correctly propagate into the steady-state
  solve — fixing a wrong-likelihood regression in estimation. A parameter
  *solved in* the steady state (``sstate_endo_param_swap``) is now seeded with a
  finite initial guess instead of leaving its NaN placeholder, and the
  parameter refresh skips NaN entries, so the solve no longer fails with
  "Objective function is returning Inf or NaN values at initial point" on
  parameter-constraint models (e.g. the canonical ``ramst``).

- **Steady-state engine replaced.** The previous (legacy) steady-state
  engine has been removed in favor of a new constraint-based engine
  built on the ``+rstst/`` package. The new engine was validated against
  the legacy one over: 55 example/OBC models (rise+solve, 0 rstst-only
  failures), order {1, 2, 3} perturbation (18/18 byte-identical),
  optimal policy including MPE (byte-identical where models parse),
  and end-to-end ``howto.m`` scripts (5/7 pass; the 2 failures are
  pre-existing bugs in the howto scripts themselves and fail
  identically under the old engine). The transient ``sstate_engine``
  option that selected between the two during validation has been
  removed; there is now only one engine. User-visible behavior should
  be byte-identical for working models; the new engine additionally
  handles occasionally-binding constraints combined with
  ``sstate_imposed`` correctly (the previous engine required a workaround).

Where the canonical changelog lives
-----------------------------------

The single source of truth for these release notes is ``CHANGELOG.md``
in the modern toolbox root, next to ``README.md``. This page is
rendered from it. New entries are added to ``CHANGELOG.md`` first, then
reflected here -- not the other way around, and never into
``rise-stable-docs``, which documents the legacy toolbox.

In MATLAB, the most recent entries are surfaced by ``rise.whatsnew``
(or its shortcut ``rise_whatsnew``):

.. code-block:: matlab

   rise.whatsnew         % print the entire CHANGELOG
   rise.whatsnew(3)      % print the 3 most recent dated/Unreleased sections
   rise_whatsnew(1)      % shortcut, same semantics

The function lives at ``+rise/whatsnew.m`` in the toolbox; the shortcut
at ``m/shortcuts/rise_whatsnew.m`` mirrors the pattern of ``rise_demo``
/ ``rise_doc``.